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19.  Abstract 


We  have  performed  detailed  theoretical  and  numerical  investigations  concerning  the 
nature  of  regions  that  arise  when  energy  is  quickly  deposited  in  a  small  domain,  i.e.,  by 
impact  or  explosion.  Such  situations  naturally  lead  to  inhomogeneous,  spatially 
dependent  probability  of  microffacturing  which  is  diminishing  with  the  distance  from 
the  center  of  the  zone  where  the  energy  was  released.  The  appearance  of  critical 
boundaries  of  fragmentation  and  of  connected  fracture  is  a  direct  consequence  of  this 
decaying  field  of  crack  probabilities. 

Within  the  fragmentation  envelope,  the  various  regions  are  distinguished  by  their 
distributions  of  fragment  sizes.  In  particular,  in  the  narrow  region  in  which 
fragmentation  just  begins  to  occur  (and  outside  which  only  connected  cracks  and 
sporadic  small  fragments  prevail),  the  fractal  form  of  the  distribution  is  found  to  be 
N~l/s'^,  where  s  is  the  fragment  volume  and  x  is  the  critical  exponent.  (For  the  two 
dimensional  case,  s  is  the  fragment  area.  The  value  of  x  varies  with  the  dimension  of 
the  space.) 

In  a  composite  target,  the  crack  probability  field  is  inherently  anisotropic. 
Generalization  of  the  classical  percolation  model  to  this  case  resulted  with  the  same 
functional  form  as  mentioned  above,  but  with  a  small  shift  in  the  value  of  x.  Size  and 
shape  distributions  of  the  fragments  in  the  2D  case  were  calculated  and  submitted  to 
the  American  partner  about  a  year  before  the  experiments  began. 

The  experimental  program,  supervised  by  Dr.  A.  H.  Mayer,  Wright  Laboratories, 
Wright  Patterson  AFB,  consisted  of  firing  1/2"  diameter  spherical  steel  projectiles  at 
carbon-fiber  epoxy-matrix  composite  plates,  made  of  32  plies.  Projectile  velocity  was 
measured  ahead  of  the  target  as  weU  as  behind  it.  Large  diameter  rigid  plastic  cylinders 
with  axes  coincident  with  the  shotline  were  lined  with  adhesive  paper  and  placed 
perpendicular  to  the  target  plate  from  both  sides.  The  fragments  ejected  from  the 
target  are  catched  on  the  adhesive  paper  and  photographed.  Each  experiment  was 
theoretically  evaluated  by  the  Israeli  group. 

Taking  the  well  recognized  view  that  near  the  ballistic  limit,  terminal  behavior  tends  to 
be  probabilistic  rather  than  deterministic,  we  defined  the  ballistic  limit  as  the  smallest 
velocity  at  which  the  front  side  fragmented  zone  (the  crater)  first  overlaps  with  the 
back  side  one  (the  spalled  zone).  That  is  to  say,  for  velocities  smaller  than  the  ballistic 
one,  these  two  fragmented  zones  are  separated  by  an  intact  layer  of  the  original  target 
material  and  therefore  the  projectile  is  not  able  to  penetrate.  On  the  other  hand,  when 
the  ballistic  limit  is  approached  from  below,  global  fragmentation  is  just  achieved.  As 
was  explained  before,  the  fractal  form  of  the  size  distribution  is  thus  expected. 

When  the  experimental  results  were  analyzed  by  the  Israeli  group  using  image 
processing  techniques,  it  turned  out  that  in  all  cases  the  fragment  size  distribution  is 
indeed  characterized  by  the  above  mentioned  functional  form.  Moreover,  the 
experimental  exponent  came  out  to  be  close  to  the  theorerical  one  (however,  subject  to 
deviations  due  to  the  finiteness  of  the  samples),  especially  for  velocities  near  the 
ballistic  limit.  Away  from  this  limit,  the  distribution  changes  (sharpens)  from  the  critical 
one  and  a  shift  in  the  experimental  exponent  is  noted  (at  higher  strain  rates,  the  relative 
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number  of  small  particles  increases).  A  microscopically  based  critical  transition 
manifestation  of  the  ballistic  limit  is  thus  established  (at  least  for  brittle  materials). 

Additional  research  efforts  during  these  years  have  been  directed  to  the  following 
subjects: 

a)  Fragment  distributions  in  a  rock  were  obtained  from  a  2D  and  3D  crack 
percolation  model  that  included  a  crack  probability  p  decaying  with  distance  from  the 
center  of  the  explosive  as  inverse  power  law.  The  mean  size  of  the  fragments  peaks 
sharply  near  the  radial  distance  r^,  at  which  the  crack  density  becomes  critical.  The 
diggable  amount  of  fragments  around  a  borehole  in  3D  exceeds  (by  several  tens  of 
percent)  the  volume  contained  inside  a  sphere  of  radius  r^.. 

b)  The  fracture  process  in  a  composite  of  short  stiff  fibers  in  a  softer  matrix 
was  simulated.  The  simulation  revealed  different  fracture  behavior  as  a  function  of 
fiber  content  and  length.  Effective  Young's  modulus,  fracture  stress  and  the  strain  at 
maximum  stress  as  a  function  of  the  volume  fraction  and  aspect  ratio  of  fibers  were 
evaluated. 

c)  Single  fiber  fragmentation  was  analyzed  using  maximum  entropy  principle. 
Theoretical  predictions  fit  quantitatively  to  experimental  ones  by  imposing  two 
constraints. 

d)  A  non-linear  transport  equation  was  derived  for  crack  population  behavior 
and  damage  induced  by  blast  in  a  medium  surrounding  a  charged  borehole.  The 
damage  is  defined  as  the  effective  volume  associated  with  the  cracks.  Intensive 
cracking  (large  damage)  appears  near  the  borehole  and  close  to  the  free  surface 
(spalling).  The  effects  of  the  initial  crack  distribution,  sample  size  and  pulse  shape  on 
the  damage  were  numerically  investigated. 
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PARTI 


ADMINISTRATIVE 
Statement  of  Work 

In  our  previous  research,  we  have  applied  percolation  theory,  Information 
theory  and  fractal  analysis  to  granular  or  porous  solids,  such  as  rocks  and  concrete,  in 
order  to  understand  basic  issues  of  damage  growth  and  fragmentation  processes.  The 
purpose  of  this  investigation  is  to  extend  our  former  work  to  advanced  materials  (such 
as  carbon-fiber  epoxy-matrix  laminates)  that  are  likewise  heterogeneous  but  also 
anisotropic. 

In  practice,  it  is  necessary  to  determine  the  critical  points  P^.  of  the  anisotropic 
system  together  with  the  universal  exponents.  Moreover,  for  fragmentation  problems 
this  is  not  enough  and  a  description  of  the  system  further  away  from  the  critical  point 
is  required  for  the  appropriate  size  (and  shape)  distributions  of  the  fragments. 

The  Soreq  group  will  take  its  share  in  the  experimental  work  in  ballistics 
carried  out  at  Wright  Laboratories,  by  proposing  discriminating  experiments  and  by 
interpreting  the  results  in  the  light  of  models  developed  at  Soreq.  The  interpretation 
will  be  aimed  at  establishing  correlation  between  the  energy  level  of  impact  and  the 
critical  densities.  In  addition,  determining  the  statistics  for  fragment  sizes  and 
comparing  them  to  the  theoretical  predictions. 

Status  of  Research 

A  mini  workshop  on  “Fracture  and  Fragmentation"  under  high  strain  rate 
loading  in  impacted  composites  was  organized  in  June  1991  at  CRREL  by  Dr.  P.K. 
Dutta  and  Dr.  A.H.  Mayer  (Wright-Patterson  A.F.B.).  Experimental  results  and 
methods  were  presented  by  the  other  participants,  while  the  Soreq  group  contributed  a 
theoretical  lecture  on  “Application  of  Percolation  Theory  to  Fragmentation  of 
Materials”. 

The  experimental  program  was  determined  as  a  result  of  this  workshop.  It  was 
planned  to  perform  a  series  of  ballistic  impact  experiments  on  plates  having  the  same 
size  and  support  with  a  projectile  having  the  same  nature  (size  and  material).  Impact 
speeds  as  well  as  rebound  velocities  (residual  ones  in  the  case  of  perforation)  are  also 
measured  in  each  experiment  in  order  to  determine  the  energy  (momentum)  absorbed 
in  the  target.  The  weight  of  the  target  is  measured  before  and  after  the  hit.  The 
weight  loss  is  a  measure  of  the  total  ejecta  mass. 

The  composite  debris  generated  by  the  impact  are  caught  at  their  point  of 
impact  with  the  duct  wall  which  is  covered  with  adhesive  paper.  The  population  of 
fragments  from  each  experiment  is  photographed  and  scanned  and  a  digitized  image  is 


formed.  Fragment  size  histograms  were  obtained  by  image  processing  techniques. 
Since  losing  fragments  (either  in  the  collection  procedure  or  due  to  low  resolution  and 
ambiguity  due  to  interference  of  shades  with  fragments’  edges)  results  in  a  biased 
fragment  distribution,  the  Soreq  group  was  actively  involved  in  this  section  of  the 
experimental  program,  as  well  as  in  the  image  processing  and  the  final  analysis. 

These  well  controlled  impact  experiments  are  tedious  and  time  consuming. 
The  first  set  of  photographs  from  each  of  the  nine  impacts  on  0°-90°  composite  plates 
arrived  at  Soreq  in  early  summer  1993.  Meanwhile,  (during  1992  and  the  first  half  of 
1993)  theoretical  models  and  tools  to  deal  with  different  aspects  of  the  problem  were 
prepared  at  Soreq.  (List  of  articles  1-5).  Especially  noteworthy  is  the  calculation  of 
fragment  size  and  shape  distribution  in  2D  anisotropic  percolation  (Part  3,  annual 
report  for  the  period  April  15,  1991-April.  14,  1992)  which  predicts  that  the  slope  of 
the  number  of  fragments  versus  fragment  size  (in  log-log  presentation)  is  r  =  2.1  (Fig. 
1,  p.  13,  above  mentioned  annual  report).  This  value  is  not  far  away  from  the 
isotropic  percolation  values,  i.e.,  x  =  2.05  (2D)  and  x  =  2.2  (3D).  Moreover,  all  these 
theoretical  predictions  came  out  very  close  to  the  experimental  slope:  x  =  2.1  (0°/90°) 
(Table  l,a,  Part  3,  this  report).  A  similar  result,  x  =  2.2,  was  accepted  from  the  set  of 
impact  experiments  on  the  0745°  composites  (Table  l,b).  More  evidence  has  been 
accumulated  since  then  (Part  4)  pointing  to  the  conclusion  that  the  Ballistic  Limit  is  a 
critical  point  in  phase  transition.  Consequently,  the  power  laws  with  universal 
exponents,  which  are  well  known  in  the  theory  of  critical  phenomena,  may  have  an 
equivalent  form  also  in  ballistics. 

The  percolation  program  was  also  extended  to  cover  fragmentation  processes 
in  three  dimensional  anisotropic  and  heterogeneous  systems.  A  new  entity,  the 
relative  excess  shattered  volume  (RESV)  was  defined.  This  quantity  is  related  to  the 
total  mass  of  material  that  is  fragmented  within  the  shattered  zone.  Using  scaling 
arguments,  we  expressed  RESV  in  terms  of  the  rate  of  decay  of  the  crack  formation 
probability  as  a  function  of  the  distance  from  the  region  of  energy  deposit. 
Simulations  in  3D  situations  show  that  RESV  increases  as  the  rate  of  decay  becomes 
milder.  In  two  dimensions,  RESV  is  practically  independent  of  the  rate  of  decay  of 
crack  formation  probabilities. 

We  have  also  treated  observationally  obtained  distributions  of  fragment  sizes 
and  showed  that  to  derive  conservation  laws  for  the  fragmentation  process  one  has  to 
compute  the  information  entropy  S  and  then  obtain  the  changes  in  S  when  the 
conserved  quantities  are  varied  by  varying  the  experimental  conditions.  The  criterion 
for  the  existence  of  conservation  laws  is  the  equality  between  the  S-derivative  and  the 
Lagrange  multiplier  in  the  Maximum  Entropy  distribution,  coupled  to  the  non¬ 
vanishing  of  the  determinant  formed  out  of  the  second  derivatives  of  S.  Application  of 
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the  method  to  ballistically  created  rock  pieces  indicates  conservation  of  the  average  of 
log  (fragment-size)  in  the  break-up  process. 

Another  topic  treated  is  the  distribution  of  broken  bonds  (as  function  of  their 
length  1)  in  a  random  2D  network,  fragmented  by  pulling.  The  obtained  distribution 
is  in  the  form  of  a  product  of  a  prior  distribution  and  an  exponential  function  (the 
constraint  factor  of  Maximum  Entropy)  whose  exponent  is  linear  in  1.  This  is  yet 
another  case  of  the  constraint  factor  being  a  simple  function  of  the  variable. 

Our  research,  in  addition,  included  approaches  to  composite  materials  other 
than  the  statistical  ones.  We  modeled  the  fiber-reinforced  materials  using  a  triangular 
lattice  of  springs  with  bond-bending  forces  between  them.  This  model  has  the 
advantage  of  being  able  to  deal  with  the  complex  situations  of  broken  fibers,  matrix 
cracks  and  fiber  matrix  debonding.  Such  situations  arise  during  the  fracture  process 
of  (continuous  or  discontinuous)  fiber-reinforced  composites.  We  first  showed  that 
the  equations  of  equilibrium  for  the  lattice  sites  reduce  to  equation  of  isotropic 
elasticity  in  the  continuum  limit.  We  then  used  the  model  to  simulate  both  the  elastic 
and  the  fracture  behaviour  of  composites.  Our  research  focused  on  the  influence  of 
fiber  volume  fraction  and  fiber  aspect  ratio  on  the  elastic  moduli  and  the  fracture 
strength  of  the  composite.  The  results  were  compared  with  other  theoretical 
approaches  and  with  experiments. 

An  attempt  to  model  the  damage  growth  and  accumulation  was  also 
accomplished.  A  non-linear  transport  equation  was  derived  for  crack  population 
behavior  and  damage  induced  by  blast  in  a  medium  surrounding  a  charged  borehole. 
The  damage  is  defined  as  the  effective  volume  associated  with  the  cracks.  Intensive 
cracking  (large  damage)  appears- near  the  borehole  and  close  to  the  free  surface 
(spalling).  The  effects  of  the  initial  crack  distribution,  sample  size  and  pulse  shape  on 
the  damage  were  numerically  investigated. 

Details  of  Work  Performed 

A  major  objective  of  the  present  research  has  been  to  generalize  the 
percolation  model  to  deal  with  the  fragmentation  process  in  highly  anisotropic 
materials  such  as  fiber  composite  laminates,  which  include  layers  of  uniaxial  high 
strength  fibers  in  a  much  weaker  matrix.  To  that  end,  we  developed  a  code  that 
assigned  different  crack  formation  probabilities  in  different  spatial  directions.  The 
code  was  applied  to  analyze  the  size  and  shape  distributions  of  the  fragments  formed 
on  lattices  of  varying  anisotropy.  A  fragment  is  defined  as  a  section  of  the  lattice 
isolated  from  the  rest  by  a  continuous  path  of  cracks.  Any  "fragment"  that  touches  the 
external  boundary  is  not  considered  in  the  distribution.  In  an  anisotropic  system,  the 
percolation  threshold  is  given  by  the  line  +  Py  =  1 .  We  have  studied  several  cases 
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on  the  threshold  line  as  well  as  a  few  cases  further  away  from  it.  A  typical  case  (Px  = 
0.7,  Py  =  0.3)  is  shown  in  Fig.  1  (page  13,  Part  3,  annual  report  for  the  period  April 
15,  1991-April  14,  1992).  The  plot  of  the  number  of  fragments  n(s)  versus  their  size  s 
on  double  logarithmic  scale  is  a  straight  line  over  at  least  3  decades  (i.e.,  between  s  = 
5  and  s  =  5000).  The  slope  that  gives  the  exponent  of  the  power  law  decay  is  about 
2.1  which  is  not  far  from  the  value  t  =  187/191  =  2.055  for  ordinary  (isotropic) 
percolation.  This  power  law  is  very  accurate  for  all  (Px,  Py)  couples  on  the  threshold 
line  within  some  range.  Outside  this  range  there  are  deviations  for  very  small 
fragments  ("dust")  as  well  as  for  very  large  ones. 

It  was  found  that  in  isotropic  systems  the  aspect  ratio  of  fragments  was  lower 
than  3.5,  while  it  could  go  as  high  as  20  in  highly  anisotropic  cases.  Photographs  of 
fragments  collected  in  the  preliminary  experiment  performed  by  P.  Dutta  at  CRJREL 
were  analyzed  using  image  processing  techniques.  The  distribution  of  fragment  sizes 
and  shapes  were  compared  with  those  obtained  from  the  simulations.  The  aspect  ratio 
distribution  from  the  experimental  fragmentation  was  found  to  have  a  width  of  about 
8,  implying  a  rather  high  degree  of  anisotropy. 

Percolation  theory  has  been  recognized  as  a  statistical  physics  model  for 
critical  phenomena.  By  applying  it  to  fracture  and  fragmentation  problems,  we  were 
able  to  show  that  fragmentation  occurs  whenever  the  density  of  the  microcracks 
reaches  its  critical  value.  That  is  to  say,  when  the  probability  p  of  creating  a 
microcrack  achieves  its  critical  value  p^.  At  and  near  this  point,  the  number  n  of 
fragments  (clusters)  having  size  s  (per  site)  behaves  as  n  oc  l/s"'^,  where  the  critical 
exponent  t  has  the  value  x  =  2.05  for  a  two  dimensional  isotropic  system  and  x  =  2.2 
for  a  three  dimensional  isotropic  one.  As  was  described  above,  we  found  the 
value  X  =  2.1  for  an  anisotropic  two  dimensional  system  (having  different  crack 
probabilities  in  the  x  and  y  orthogonal  directions).  These  are  the  statistical  physics 
predicted  results  for  the  fragment  size  distribution  at  the  critical  point,  which  is 
identified  as  the  ballistic  limit  (see  Part  4). 

The  experimental  program  performed  at  W.P.A.F.B.  is  briefly  described  in 
Part  4.  Additional  details  relating  to  the  image  analysis  procedure  are  given  in  Part  3. 
Consider  the  trend  of  the  damage  and  the  energy  transfer  ratio,  defined  here  to  be  the 
ratio  of  imparted  to  initial  projectile  energies,  as  the  projectile  velocity  (or  energy)  is 
increased.  At  very  low  velocities,  the  projectile  experiences  a  perfect,  loss-free 
rebound  from  the  panel  and  the  transfer  ratio  remains  at  zero.  As  the  velocity  is 
increased  beyond  the  damage  threshold,  a  crater  is  produced  in  the  panel,  causing  the 
transfer  ratio  to  increase.  With  higher  impact  velocities,  incipient  spall  first  appears 
and  the  transfer  ratio  further  increases.  Below  the  ballistic  limit,  it  is  expected  that  the 
spalled  fragmented  zone  at  the  backside  of  the  sample  and  the  shattered  and 
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penetrated  zone  in  the  front  side  (the  crater)  are  separated  by  an  intact  or  delaminated 
layer  of  the  material  which  is  strong  enough  to  prevent  the  projectile  from  perforating 
the  target.  In  the  limit  (in  fact,  the  ballistic  limit)  this  intact  layer  becomes  very  thin 
and  the  two  damage  zones  will  touch  each  other.  That  is  to  say,  the  sample  will  then 
become  disintegrated  throughout  its  thickness. 

Our  experimental  investigation  provides  quantitative  evidence  that  the  ballistic 
limit  represents  a  critical  phenomenon.  This  evidence  comes  from  the  size 
distribution  of  the  ejecta  produced  near  the  ballistic  limit.  The  distribution  density 
function  for  fragments  of  area  s  was  found  to  decrease  with  s  raised  to  a  negative 
power.  Results  for  several  impact  velocities  for  0°-90°  fiber  orientation  are  presented 
in  Table  l,a.  Similar  results  for  0°-45°  are  given  in  Table  l,b  (Part  3).  The 
experimental  exponents  t  =  2.1  (0°-90°)  and  t  =  2.2  (0°-45°)  are  very  close  to  the 
values  found  theoretically  for  fragment  sizes. 

Impact  loading  or  other  localized  energy  deposition  mechanisms  (i.e. 
explosion  and  laser  irradiation)  lead  to  spatially  dependent,  inhomogeneous 
probability  of  microfracturing.  A  three  dimensional  percolation  simulation  of 
fragmentation  has  been  performed  for  the  first  time  (at  this  stage  only  for  isotropic 
microfracturing.)  The  characteristic  (average)  fragment  size  as  a  function  of  the 
radius  from  the  point  at  which  the  energy  was  released  was  obtained  in  two  and  three 
dimensions.  Close  to  the  impact  point  fine  fragments  are  obtained,  but  their  size 
increased  up  to  a  maximum  just  short  of  the  distance  from  the  center,  which  is 
defined  as  the  position  at  which  the  crack  formation  probability  achieves  the  critical 
value  pc^I  for  fragmentation.  The  distance  R^  is  also  meaningful  for  estimating 
quantitatively  the  total  amount  of  shattered  material.  In  a  naive  picture,  the  volume 
(surface  in  2D)  of  the  broken  material  is  just  equal  to  4/37r(RcI^)2  in  2D).  Deviations 
from  these  values  have  been  calculated  for  various  energy  deposition  configurations 
and  were  found  to  be  the  more  in  access  of  the  naive  values  the  shallower  is  the 
decrease  of  the  microcracing  probability  with  the  distance  from  the  center.  (Part  2). 

Fragmentation  is  a  complex  phenomenon.  Its  outcome,  the  fragment 
distribution,  is  expected  to  reflect  this  complexity.  The  Maximum  Entropy  (ME) 
approach  has  been  used  in  the  past  to  impose  regularities  on  the  distribution.  In  this 
work  we  set  out  to  uncover  regularities  from  the  data.  The  key  to  this  effort  are  the 
derivatives  of  the  entropy  with  respect  to  the  invariants.  The  procedure  leads  to  a 
"Thermodynamics”  of  the  fragment  -  (or  any  other)  data.  Application  has  been  made 
to  relatively  simple  case  of  ballistic  fragmentation  of  rock  targets,  where  only  one 
invariant  <logs>  is  in  evidence.  For  future  uses  of  the  procedure,  when  more  than  one 
invariant  is  expected,  further  development  of  the  method  is  needed  to  make  it 
practical.  This  appears  possible.  (Part  5). 
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A  priori  -  constraint  function  decomposition  of  distributions  has  been  also 
carried  out  in  this  work  for  the  broken  bonds  in  a  disordered  network  at  the  moment  of 
its  falling  to  pieces.  The  aim  was  to  see  whether  the  constraint  function  is  simple,  as 
anticipated  in  the  ME  treatment.  This  was  found.  Additional  interest  is  attached  to 
possible  changes  in  the  distribution  as  the  fragmentation  threshold  is  approached. 
(Part  5). 

The  transport  equation  approach  to  crack  population  growth  (Part  6)  is  only 
one  of  the  possible  treatments  of  a  complicated  stochastic  problem:  to  us  it  appears 
neither  less  nor  more  justified  than  any  other,  and  for  a  space  and  time-dependent 
description  it  may  be  uniquely  workable.  There  are,  admittedly,  a  great  number  of 
free  parameters  in  the  approach,  whose  values  are  not  easily  fixed  in  any  concrete 
situation.  Our  aim  was  to  elucidate  the  wide  range  of  phenomena  that  can  occur  with 
some  (not  unphysical)  choice  of  the  parameters.  Yet  more  parameters  can  be  varied 
(as  the  crack  velocity  or  introducing  inhomogenity  in  the  initial  crack  distribution). 

Although  the  results  presented  in  his  work  are  related  to  the  specific  samples 
defined  in  Part  6,  it  is  possible  to  reach  some  conclusions,  which  seem  to  present  the 
general  behavior  of  damage  in  systems. 

(a)  In  a  given  sample,  a  pressure  wave  creates  a  spall  zone,  if  the  ratio  is 

larger  than  some  threshold  value. 

(b)  For  a  given  pressure  wave  in  a  given  medium,  the  damage  peak  in  the  spall 
zone  increases  for  Rq  <  R*o  and  decreases  for  Rq  <  R*o-  For  sufficiently 
large  size  of  the  sample  Rq,  the  spall  does  not  exist.  (In  the  case  of  Sec.  4,  the 
threshold  is  R*o  =  4.4). 

(c)  For  a  given  pressure  wave  in  a  given  medium,  the  damage  peak  and  the  width 
of  the  inner  damage  zone  increases  with  the  sample  radius  Rq. 

(d)  For  a  given  pressure  wave  and  sample  size  Rq,  if  a  >  a*  is  a  sufficiently  small 
number,  the  damage  peak  and  the  width  of  the  inner  damage  zone  and  the  spall 
region  are  independent  of  the  a  value.  In  the  limit,  as  a  ^  0,  the  peak  of  the 
damage  near  the  center  reaches  a  minimal  value,  and  the  peak  of  the  spall 
attains  a  maximal  one. 

(e)  For  a  given  pressure  wave  and  sample  size  Rq,  if  a  >  a*,  where  a*  is 
sufficiently  large  number,  the  spall  zone  does  not  exist. 

(f)  A  pressure  wave  and  a  sample  of  radius  Rq  are  given.  For  5  >  5*,  is  a 
sufficiently  large  number,  the  spall  zone  does  not  exist  and  the  damage  peak  of 
the  inner  damage  zone  increases  with  the  5  value.  In  the  limit,  5  — >•  oo,  the 
damage  achieves  a  maximal  value. 


7 


Participating  Personnel 

Robert  Englman  (Prof.);  B.Sc.,  A.R.C.S.  (in  Mathematics)  in  1954  and  Ph.D. 
(Theoretical  Physics)  in  1957  at  London.  Thesis  title  "Extension  of  the  theory  of  the 
mean-free-path  effects  in  metals". 

Zeev  Jaeger  (Dr.);  M.Sc.  in  Chemistry,  Physics  as  major  subjects  and  Mathematics 
and  Applied  Mathematics  as  subsidiary  subjects  in  1967  at  the  Hebrew  University  in 
Jerusalem,  Ph.d.  in  Physics  in  1973  at  the  Hebrew  University,  Jerusalem  with  the 
thesis  title.  "Charge-transfer  states  and  their  physical  application  in  crystals  of 
transition  metals". 

Raphael  Ruppin  (Dr.);  M.Sc.  in  Physics  at  the  Hebrew  University  of  Jerusalem  in 
1964  and  Ph.D.  in  Theoretical  Solid  State  Physics  from  the  Hebrew  University  in 
1969  with  the  thesis  title  "Size  dependent  elementary  excitations  in  solids". 

Miriam  Lemanska  (Dr.);  M.Sc.  in  Mathematics  in  1952  at  the  Lodz  University  in 
Mathematics  and  Ph.D.  in  Nuclear  Sciences  in  1967  at  the  Technion,  Haifa  in  the 
subject  "Perturbation  theory  in  fast  reactors". 

Michael  Murat  (Dr.);  B.Sc.  in  Chemistry  in  1977  from  the  Technion,  Haifa,  then 
M.Sc.  in  Physical  Chemistry  in  1978  from  the  Technion,  Haifa  with  a  thesis  on  "The 
triplet  state  of  the  biexciton"  and  Ph.D.  in  Physics  in  1987  at  the  Tel  Aviv  University 
with  a  thesis  entitled  "The  structure  and  physical  properties  of  porous  materials". 

Micha  Anholt  (Mr.);  B.Sc.  in  Physics  (1989)  from  Tel  Aviv  University.  Studying 
towards  M.Sc.  in  Physics  of  Condensed  Matter. 


Interactions 


Prof.  T.J.  Ahrens,  Dept  of  Geology  and  Astrophysics,  California  Institute  of 
Technology,  Pasadena.  CA,  Discussions  with  Z.  Jaeger  July  1991,  June  1992. 

Prof  H.D.  Wagner,  Department  of  Materials  and  Interf  The  Weizmann  Institute  of 
Science,  Rehovot  76100,  Israel.  Discussion  with  M.  Murat  throughout  the 
year. 

Dr.  E.  Schulz  &  Dr.  A.  Hampe,  Federal  Institute  for  Materials  Research  and  Testing 
(B.A.M.)  D-1000,  Berlin  45,  Germany.  Correspondence  during  the  year  and 
discussions  with  R.  Englman  and  M.  Murat  in  Israel,  June  1992. 

Prof  S.L.  Phoenix,  Dept,  of  Theoretical  and  Applied  Mechanics  Cornell  Univ.,  Ithaca, 
N.Y.  Discussions  with  Z.  Jaeger  at  Cornell. 

Dr.  J.  Karger-Kocsis,  Institute  for  composite  materials,  D-6750  Kaiserslautern, 
Germany.  Correspondence  during  the  year  and  discussions  with  R.  Englman 
and  M.  Murat  in  Israel,  June  1992. 

Prof  S.M.  Wiederhom,  National  Institute  of  Standards  and  Technology.  Several 
discussions  with  Z.  Jaeger  1992-1993. 

Dr.  F.  Allahdadi,  Space  Kinetic  &  Impact  Debris  Group,  Space  Technical  Center 
Philips  Laboratories,  Albuquerque,  N.M.  Presentation  by  Z.  Jaeger  and 
discussion,  1993. 

Dr.  D.E.  Grady,  Experimental  Impact  Physics  Div.  1433,  Sandia  N.L.,  Albuquerque, 
N.M.  Correspondence  and  discussions  1992-1995. 

James  F.  O’Bryon,  Deputy  Director,  Operational  Test  and  Evaluation,  Live  Fire 
Testing,  Office  of  the  Secretary  of  Defense  DDOT  &  E/LFT,  The  Pentagon. 
A  few  discussions  on  damage  in  composites  and  vulnerability  of  flying 
structures  made  of  composites,  Jerusalem,  Israel,  May  1955. 


9 


List  of  Articles 


1.  R.  Englman  and  Z.  Jaeger,  “Percolation  Theory  of  Fragmentation”,  an  invited 
lecture  presented  at  the  International  Symposium  on  Intense  Dynamic  Loading 
and  its  Effects,  Chengdu,  China,  June  9-12  (1992)  pp.  425-430. 

2.  M.  Murat,  M.  Anholt  and  H.D.  Wagner,  ‘Tracture  Behaviour  of  Short-Fiber 
Reinforced  Materials”.  J.  Mater.  Res.  7, 3120-3131  (Nov.  1992). 

3.  R.  Englman,  “Information-Entropy  Distribution  for  Sequential  Fragmentation 
in  Single  Fibers”,  Int.  Conf.  on  “Microphenomena  in  Advanced  Composites”, 
G.  Marom  and  H.D.  Wagner  (Edits)  Herzlia,  Israel,  June  28- July  1,  (1992). 

4.  M.  Anholt,  R.  Englman  and  Z.  Jaeger,  “Percolation  Approach  to  Locally 
Caused  Fragmentation”,  Centre  de  Physique  des  Houches  Workshop  on 
Fragmentation  Phenomena,  12-17  April  1993,  Editor:  D.  Beysens. 

5.  Z.  Jaeger,  M.  Anholt  and  A.H.  Mayer”,  Impact  of  Fiber  Composite  Laminate 
Plates:  A  Percolation  View  of  Perforation  and  Spallation”  in  “Dynamic 
Response  of  Composite  Structures”,  D.  Hui  and  G.  Anderson  (Edits),  New 
Orleans,  LA,  August  30-Sept.  1, 1993. 

6.  A.H.  Mayer,  Z.  Jaeger  and  R.  Englman,  “The  Ballistic  Limit  of  Laminated 
Carbon  Fiber  Epoxy  Plates  as  a  Critical  Point  in  Phase  Transition”,  a  post 
deadline  paper  (TB60)  in  the  15th  International  Symposium  on  Ballistics, 
Jerusalem,  Israel,  21-24  May  1995.  M.  Mayseless  and  S.R.  Bodner  (Edits.). 

7.  R.  Englman,  “Maximum  Entropy  Principles  in  Fragmentation  Data”  in  “High 
Pressure  Shock  Compression  in  Solids  II”.  Editors:  L.  Davison,  DE.  Grady 
and  M.  Shahimpoor,  Springer-Verlag,  N.Y.  (1994). 

8.  M.  Lemanska,  R.  Englman  and  Z.  Jaeger,  “Transport  Treatment  of  Crack 
Population  in  Finite  Medium  Part  11:  Damage  Behavior  in  Different  ZIbnes”. 
Submitted  for  publication. 


Patents  and  Inventions 


None 


r 


10 


Statement  to  Program  Manager 
Achievements 

The  first  summer  of  research  began  with  presentations  of  earlier  achievements  which 
were  presented  by  the  participants  of  the  mini  workshop  on  composites'  impact 
problems.  As  a  result  of  the  preliminary  impact  experiment,  the  lectures  and  the 
discussions  during  the  workshop,  the  directions  of  research  during  the  year  were 
agreed  upon. 

During  the  second  year  percolation  model  was  generalized  to  anisotropic  two 
dimensional  case  as  well  as  to  the  nonhomogeneous  and  anisotropic  three  dimensional 
one.  Shape  and  size  distributions  of  fragments  in  the  2D  case  were  calculated  and 
compared  at  first  with  data  deduced  by  image  processing  techniques  from  results 
achieved  by  the  preliminary  experiment.  A  series  of  very  well  controlled  impact 
experiments  with  different  velocities  on  0790°  composite  plates  were  performed  at 
W.P.A.F.B.  and  analyzed  several  months  after  the  annual  report  for  the  second  year, 
(including  the  theoretical  percolative  values)  was  submitted.  It  was  found  that  the 
distribution  density  function  for  fragments  of  area  s  decreases  with  s  raised  to  a 
negative  power.  Moreover,  the  value  of  the  experimental  negative  power  at  the 
ballistic  limit  was  indistinguishable  from  the  theoretically  predicted  one.  Later  on  a 
second  series  of  similar  experiments  was  done  on  0°/45°  composite  plates  and  the 
analysis  gave  again  very  similar  results.  The  evidence  accumulated  show  that  we 
have  succeeded  in  uncovering  the  physical  meaning  of  the  perforation  threshold 
(the  Ballistic  Limit):  it  represents  a  critical  phenomenon.  Moreover,  Percolation 
theory  could  provide  an  adequate  description  of  the  fragmentation  process.  The 
identification  of  the  ballistic  limit  velocity  with  a  critical  point  lays  open  the  field 
of  ballistics  to  the  theoretical  apparatus  of  critical  processes  and  its  scaling  laws.  In 
addition,  it  has  also  been  shown  that  much  insight  into  the  complicated  process  of 
fragmentation  can  be  achieved  by  using  the  Maximum  Entropy  Principle. 

Future  Research  Envisaged 

The  fragmentation  process  is  clearly  related  to  the  damage  accumulated  in  the 
impacted  composite  at  and  above  the  Ballistic  Limit.  It  is  expected  that  scaling  laws 
which  hold  at  the  critical  point  of  the  phase  transition  will  also  hold  around  it.  New 
hitherto  unknown  scaling  laws  for  damage  would  be  theoretically  identified  and 
experimentally  confirmed.  Such  relations  are  relevant  to  vulnerability  considerations 
of  flying  systems. 

It  is  suggested  to  repeat  the  experiments  with  velocities  near  and  around  the 
ballistic  limit.  It  is  recommended  to  use  better  techniques  of  photographing  and 
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higher  resolution  equipment  for  image  processing  in  order  to  receive  more  reliable 
experimental  fragment  size  distributions  (fragment  size  should  span  two  or  three 
orders  of  magnitude).  It  is  suggested  that  a  series  of  well  controlled  impact 
experiments  would  be  repeated  with  different  target  thicknesses,  and  different  sizes 
and  shapes  of  projectiles.  (Control  on  the  strain  rate  could  be  achieved  by  changing 
the  sphere  to  a  long  cylinder).  Larger  scale  experiments  should  also  be  accomplished. 
Ejecta  at  the  front  side  has  to  be  measured  and  analyzed  along  the  same  lines  that 
ejecta  from  the  backside  was  studied.  It  is  also  recommended  to  trace  back  spacial 
effects  in  damage  and  fragmentation  by  identifying  the  places  from  where  the 
fragments  and  their  velocities  came  from. 


PART  2 


PERCOLATION  APPROACH  TO 
LOCALLY  CAUSED  FRAGMENTATION 


Micha  Anholt,  Robert  Englman  and  Zeev  Jaeger 


Department  of  Physics  and  Applied  Mathematics 
Soreq  Nuclear  Research  Center 
Yavne  70600,  Israel 


ABSTRACT:  The  paper  starts  with  a  survey  of  the  variety  of  fragmentations, 
with  emphasis  on  the  factors  which  determine  whether  the  size  distribution  will  be 
monotonic,  decrecLsing  with  size  or  peaked.  We  consider  cases  of  loading  that  are 
local  and  uniform,  fast  and  slow;  confined  or  locked-in  situations  and  excess-energy 
redistributions  that  are  channel-specific  or  statistical.  Fragment  distributions  in 
a  rock  were  obtained  from  a  2D  and  3D  crack  coalescence  (percolation)  model, 
that  included  a  crack  formation  probability  p  decaying  with  the  distance  from  an 
explosive  bore-hole  as  an  inverse  power  law.  The  mean  size  of  the  fragments  peaks 
sharply  near  that  radial  distance  Vc  at  which  the  crack  density  is  becomes  critical. 
In  this  neighbourhood  the  spread  and  the  height  of  the  mean  fragment-size  peak 
follow  scaling-laws  appropriate  to  space- varying  p.  The  digable  amount  of  fragments 
around  a  borehole  in  three  dimensions  exceeds  (by  several  tens  of  percents)  the 
volume  contained  inside  a  sphere  of  radius  Tc  . 
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Introductory  Overview 

The  subject  of  fragmentation  is  exceedingly  ramified  and  calls  for  classification. 
This  is  attempted  here  from  several  viewpoints. 

First  we  make  some  rather  obvious  classifications  related  to  the  mode  of  load¬ 
ing  or  of  energy  input.  This  can  be  done  in  a  manner  affecting  the  medium  uni¬ 
formly  (such  as  occurs  when  one  bends  a  bar)  or  in  a  local  manner  (which  is  the 
case  treated  in  the  following  sections  in  this  paper).  In  the  former  the  stress  in 
the  medium  essentially  scans  the  medium  for  its  weakest  points  (and  is  strain  con¬ 
trolled)  ,  while  in  the  latter  fracture  may  occur  at  the  point  where  stress  is  applied 
and  propagates  from  that  point.  The  speed  of  the  loading  has  essentially  two  con¬ 
sequences  on  fragmentation:  It  affects  characteristic  particle  sizes  [1,2,3]  ,  in  the 
sense  that  high  strain  rates  produce  small-size  fragments  and  tend  to  act  uniformly 
through  the  medium.  The  number-distribution  of  fragment-sizes  found  experimen¬ 
tally,  or  predicted  theoretically  often  peaks  at  a  finite,  nonzero  size,  the  weight- 
weighted  distribution  (important  in  technology)  almost  invariably  does  so.  This 
statement  holds  for  one-time,  sudden  loading;  sequential  breakage  favours  small 
sized  fragments  whose  distribution  has  some  inverse  power  law  form  (fractality  [4]). 
An  essentially  theoretical  postulate  affects  directly  the  predicted  distribution.  If 
one  assumes  that  the  loading  energy  is  a  priori  equally  likely  distributed  among  all 
bonds  between  particles,  then  a  peaked  distribution  is  expected;  on  the  other  hand, 
if  all  fragments  are  estimated  to  be  a  priori  on  equal  footing,  then  a  monotonic 
distribution  (decreasing  at  larger  sizes)  results.  Both  cases  are  quoted  in  a  review 
[3];  however,  it  is  not  clear  which  choice  fits  which  factual  circumstance.  Compe¬ 
tition  for  the  available  energy  by  different  modes  has  remarkable  consequences  in 
molecular  physics  [e,g.y  photofragmentation  and  predissociation  [5])  and  presum¬ 
ably  in  nuclear  multi-fragmentation.  When  excitation  is  brought  about  locally  the 
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spread  of  excitation  energy  in  time  and  space  is  a  subject  that  has  been  researched 
especially  since  the  70’s.  In  one  extreme  there  is  the  possibility  of  a  fast  redistri¬ 
bution  of  the  excitation-energy  among  all  modes  democratically  (this  is  the  case  of 
thermalization  in  the  ‘statistical’,  high  density-of-state  limit )  which  tends  to  result 
in  the  simultaneous  production  of  many  molecular  fragments.  Alternatively,  one 
has  situations  (for  harmonically  bonded  particles)  that  the  mode  with  the  highest 
frequency  takes  up  most  of  the  energy  and  breaks  preferentially.  At  the  other  ex¬ 
treme  the  energy  becomes  channelled  to  a  special,  localized  (highly  non-harmonic) 
mode,  which  is  the  one  that  breaks  first,  while  the  other  modes  are  only  moderately 
excited.  A  further  classification  depends  on  whether  fragmentation  occurs  when 
the  matter  is  in  a  confined,  locked-in  situation  or  it  is  free.  (Rock  blasting  and 
spacecraft  breaking  up  are  representative  instances  for  the  two  cases.)  There  are 
indications  that  the  latter  gets  fragmented  by  a  percolation  process.  This  means 
that  cracks  are  formed  fairly  rapidly  and  uniformly  in  the  system  and  when  these 
coalesce  the  fragments  fly  apart  without  further  interaction  with  each  other.  The 
coalescence  is  a  percolation  process.  The  fragment  size  distribution  of  space  debris 
appears  to  arise  under  such  circumstances  [3]. 

Background 

Rock  blasting  for  the  purpose  of  mining  is  carried  out  by  exploding  point-like  or 
cylindrical  charges  in  bore-holes.  As  a  consequence,  cracks  are  created  whose  density 
near  the  borehole  is  large  and  there  the  rock  breaks  into  small  fragments  (the 
crushed  zone).  As  one  moves  away  from  the  centre  the  crack-density  decreases  and 
larger  fragments  start  to  dominate  (the  fragmented  zone) ,  until  one  reaches  region 
where  the  ‘infinite-fragment’  contains  mainly  isolated  cracks  (the  cracked  region), 
ending  in  intact,  uncracked  recesses. 
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The  Simulational  Scheme 

In  this  work  the  statistics  of  the  fragment  in  obtained  and  analyzed  in  the  framework 
of  a  percolation  model,  that  extends  to  three  dimensions  the  earlier  computations 
of  Murat  et  al.  [6].  The  rock  is  modelled  by  a  cubic  lattice(or  by  a  square  lattice 
in  2D)  and  planar  crack  plaquettes  are  formed  on  the  faces  of  the  cells  of  linear 
size  a,  (representing  the  granular  composition  of  the  rock)  with  a  probability  which 
decreases  with  the  distance  r  from  the  centre  according  to 

pW  =  (’■o/r)“  ;  I’'  >  ’■o)-  (1) 

ro  is  the  borehole  radius,  a.  =  0. 2-0.4. 

We  expect  from  percolation  theory  for  infinite  homeogeneous  structures  [7] 
that  at  distances  below  the  value  where  p{r)  =  pd  (the  standard  critical  percolation 
probability;  for  numerical  values  of  this  and  other  critical  parameters  in  2D  and  3D, 
see  table  I)  there  will  be  a  connected  cluster  of  cracks  across  which  (for  example) 
fiuid  can  percolate  throughout  the  inside  volume.  However,  as  shown  in  Ref.  [8], 
fragmentation  will  not  occur  inside  all  this  volume,  but  only  within  a  smaller  volume 
defined  by  the  critical  fragmentation  radius  rc  at  which 

p{rc)=Pc2  (>Pcl)- 

see  Figure  1.  At  distances  less  than  rc  fragmentation  will  be  intensive  and  the  mean 
fragment  size  will  be  small,  approaching  the  minimal  size  a  (the  ‘grain’  dimension). 
For  distances  beyond  one  encounters  the  infinite  fragment  plus  some  small  frag¬ 
ments  embedded  in  it.  Due  to  considerations  that  are  both  practical  and  theoretical, 
we  exclude  henceforth  the  infinite  fragment  from  summation  over  fragments.  Then 
at  r  >  Tc  the  mean  fragment  size  will  decreaise,  once  again  approaching  a,  since  at 
low  crack  densities  fragments  formed  out  of  more  than  the  minimal  crack  number 
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(4  in  2D  ,6  in  3D)  are  extremely  unlikely.  Mean  fragment  size  Sav  a-s  function  of  the 
radius  are  shown  in  Figures  2,3.  Computations  were  performed  on  lattices  of  sizes 
up  to  (I50a)^  -  in  3D  and  (2000a) ^  in  2D,  sufficiently  large  to  neutralize  any  edge  or 
finite  size  effects.  There  were  between  50  (in  3D)  and  2000  (in  2D)  realizations  for 
each  set  of  parameters.  The  exponent  a  in  (l)  was  taken  to  vary  between  0.2~0.4,2is 
already  remarked;  the  borehole-variation  is  discussed  below. 

Scaling 

Since  the  total  lattice  size  is  large  enough  not  to  affect  the  results  (els  is  apparent  in 
Figure  1),  there  remain  only  two  length  variables:  the  borehole  radius  tq  and  the 
cell  size  a.  In  making  the  model  apply  to  an  actual  blast,  we  assign  to  the  borehole 
radius  a  value  that  corresponds  to  the  actual  situation  ,  Tc  is  then  obtained  (at 
least  roughly)  from  the  outcome  and  will  be,  if  the  model  is  anything  like  right, 
of  the  order  of  tq,  independent  of  a.  However  the  average  size  Sav  Figures  2,3 
and  the  width  of  the  distribution  are  given  in  units  of  a  and  the  effect  upon  these 
of  increasing  tq  is  equivalent  to  a  corresponding  decrease  in  a.  This  circumstance 
is  useful  since  we  have  varied  the  borehole  radius  tq  and  computed  changes  in  the 
location  of  the  maximum  of  5at;  and  in  the  width  w  of  the  peak  in  the  Sav  versus  r 
curve(Figure  4).  It  is  apparent  that  the  location  of  the  5av-niaximum  scales  linearly 
with  the  borehole  size,  which  is  reasonable  since  its  value  is  close  to  (though  not 
quite)  Tc  and  this  scales  linearly  with  rg.  In  the  same  figure  the  width  shows  a  milder 
variation  with  rg  which  is  linear  with  rg  for  large  boreholes  and  something  like  Tq*^ 
for  small  boreholes.  Recalling  our  previous  remark  regarding  reciprocity  between 
rg  and  a  we  can  formalize  this  bs  w fa  oc  a“^*^  for  large  -a,  or  equivalently  w  a  a^*^. 
To  discuss  quantitatively  behaviors  of  Sav  and  of  w  in  the  critical  regime,  namely 
around  r  =  rc,  we  recall  that  in  a  lattice  with  uniform  occupation  probability  both 
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the  maximum  fragment  size  and  the  width  will  be  of  the  order  of  the  size  L  of 
the  lattice.  For  a  non-uniform  probability  the  spatial  variation  of  p{r)  around  Tq 
will  provide  the  appropriate  size  scale  and  it  has  been  shown  [6,9]  that  the  correct 


scaling  is  with 


Leff  «  {'rc/aPc2) 


I'/il+u) 


We  shall  now  relate  the  average  size  the  fragments  near  Tc  to  the  maximum 
fragment  size  Smax  utilizing  the  result  [3]:  Ug  oc  s~^,  where  Ug  is  the  average 
number  per  site  of  fragments  having  size  5.  (We  assume  the  statistical  isotropy  of 
fragments  though  the  model  does  not  warrant  this.)  Then 


Smaxja'^  ,  rSmax/a'^ 

dsUgS^  I  J  dsTlgS 


.T  -  2. 


Using  the  values  in  Table  I  this  leads  to  the  following  scaling  between  5at;(near 


r  =  Tc)  and  a 

Sav  oc  in  2D 

Sav  OC  in  3D 

Let  it  be  remarked  that  away  from  Tq  on  both  sides  Sav  varies  with  a^.  Our 
numerical  results  confirm  this.  To  calculate  the  width  w  of  the  Sav  peak  around 

we  remark  that  this  is  related  to  the  fractal  dimension  by 

^  _  / Smax\\l D 

a~^  ad  ’  (4) 

From  the  first  line  we  obtain  w  =  [Smax)^^^ and  from  the  second  line 
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w  oc  3D 

which  may  be  compared  with  the  exponent  0.5  estimated  from  computation..  Nu¬ 
merical  comparison  (Table  II  )  shows  that  there  is  a  constant  factor  of  4  between 
the  two  members  of  Eq.  (5)  in  both  2D  and  3D.  (Presumably  this  means  that  the 
definition  of  the  width  is  not  clear  cut)  Figure  5  exhibits  the  probability  poo  l^hat  a 
site  belongs  to  the  infinite  fragment  as  fuction  of  the  distance  from  the  centre.  The 
objective  is  to  quantify  what  part  of  the  rock  can  be  mined  or  digged  out,  with  the 
understanding  that  the  infinite  fragment  is  not  digable.  In  the  uppper  curve  the 
isolated  or  embedded  fragments  are  counted  as  part  of  the  infinite  fragment  and 
the  approach  to  unity  is  faster  than  when  the  embedded  fragments  are  presumed 
digable  (the  lower  curve).  Estimation  of  digability  is  the  subject  of  the  following 
section. 


The  Excess  Digable  Fraction  (EDF). 

Digable  broken  pieces  appear  beyond  because  fragments  have  finite  size  on  the 
scale  of  the  basic  crack  size  a,  and  are  not  infinitesimal,  point  like.  Likewise, 
unbroken  pieces  belonging  to  the  infinite  fragment  appear  inside  Tc.  Both  effects 
represent  fluctuations  from  an  ‘ideal’  situation  in  which  the  radius  r  —  Tc  defines 
a  sharp  boundary  between  fragmented  and  unbroken  (though,  possibly,  cracked) 
regions.  Both  small  and  large  sized  fluctuations  exist  and  are  apparent  in  Figure  1  ; 
however  statistically  the  most  important  ones  are  the  large  sized  ones.  The  digable 
volumes  as  function  of  fragments-volume  are  shown  in  Figures  6-8  for  3D  simulations 
with  the  charge  having  either  spherical  or  cylindrical  geometries.  Parameters  are 
marked  in  the  figures.  The  horizontal  line  represents  the  ideal  situation  where  r  = 
(written  as  Rq  in  the  drawing)  is  a  sharp  boundary.  We  now  define  the  excess  digable 
fraction  (EDF)  as  the  relative  difference  between  the  computed  digable  volume  and 
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the  ideal  one.  The  excess  fraction  is  also  the  difference  between  averaged  volume  of 
broken  pieces  (including  those  that  extrude  beyond  fc)  and  the  averaged  volume  of 
the  unfragmented  pieces  (including  the  intrusions),  all  this  divided  by  the  volume 
contained  within  Figure  9  represents  the  results  of  our  computation  in  2D  and 
in  3D  (with  two  types  of  charges).  Now  in  two  dimensions  on  a  square,  where  the 
opened  and  closed  plaquettes  are  self-dual  (with  Pc2  =  0.5  for  either),  fragments 
and  unbroken  pieces  behave  identically  and  one  expects  that  to  within  random 
statistical  error  there  will  be  no  net  excess  fraction.  In  three  dimensions  EDF  will 
have  the  form 

EDF  (X  (L,/ cc  pl3-/’‘^/(l+‘^)l/“_  (5) 

We  have  used  here  and  have  inserted  in  the  exponent  the  parameter 

p  whose  value  is  between  1  and  3,  due  to  our  uncertainty  cls  to  the  shape  of  the 
extracted  volume  (whether  it  is  similar  to  a  sphere  ,  yo  =  3,  or  to  a  flat  piece  lying 
on  the  hull  at  p  =  1).  In  any  case,  due  to  both  factors  in  EDF,  this  decreases  with 
increasing  value  of  the  exponent.  In  our  simulations  a  positive  EDF  is  found  in  3D 
and  a  power  law  dependence  in  the  form  EDF  =  0.05a“^'^  .  Contrasting  this  are 
the  two  dimensional  results  in  Figure  9  which  show  that  EDF  is  negative  and  so 
small  as  to  be  zero  within  the  statistical  errors  of  our  data;  The  dependence  on 
a  is  also  negligibly  small.  These  agree  with  the  analysis  of  the  2D  situation  given 
above. 
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Conclusion 

Our  model  of  fragmentation  based  on  the  coalescence  of  randomly  created  cracks 
needs  a  good  deal  of  experimental  confirmation,  especially  when  applied  to  a  con¬ 
fined  geometry  (as  in  solid  rock).  However,  in  fast  processes  which  characterize 
blasting  there  is  some  justification  for  the  model,  including  the  neglect  of  correla¬ 
tions  in  crack  development.  Though  our  results  are  based  on  a  special  form  [Eq.(l)] 
of  crack  generation  probability,  there  are  several  generally  valid  features  which  could 
could  be  useful  in  the  design  of  post-blast  mining  operation  (quarrying).  Examples 
are  the  rather  sharp  maxima  in  fragment  sizes  near  Tc,  the  advantages  we  have  found 
in  having  sharply  decaying  probabilities  (a  large) relative  to  fiat  ones  (a  small),  and 
finally  our  finding  that  there  are  different  digable  volumes  in  different  dimensions. 
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Table  I 

Values  of  parameters  used  in  our  computation(symbols  explained  in  text).  Critical 
exponents  and  concentrations  rely  on  Ref  [7]. 


2D 

3D 

d 

2 

3 

D 

1.9 

2.5 

Z-r 

0.9 

0.8 

V 

1.33 

0.9 

Pci 

0.5 

0.249 

Pci 

0.5 

0.751 
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Table  II 

Maximum  value  of  average  size  Sav  near  Tc  and  width  w  of  the  distribution  of  average 
sizes  as  function  of  distance  from  the  borehole  obtained  from  our  computation 
(expressed  in  units  of  the  cell  dimension  a). 


2D 

3D 

960 

38 

4(Sav/a^?l^ 

148 

17.1 

[wlaY 

150 

17 

*  From  computation 


12 


24 


FIGURE  CAPTIONS 


Figure  1:  One  realization  on  a  square  lattice.  Different  clusters  are  marked  by 
different  grey  shades.  The  black  lines  represent  the  cracks  opened  cracks. 

Figure  2:  Sav  ^  function  of  r  for  a  square  lattice  of  linear  size  2000a,  borehole 
radius  100a,  a  =  0.35.  Results  were  averaged  over  50  realizations. 

Figure  3:  Sav  as  function  of  r  for  a  cubic  lattice  of  linear  size  150a,  borehole 
radius  10a,  a  =  0.2.  Results  were  averaged  over  200  realizations. 

Figure  4:  Location  of  maximum  of  Sav  and  width  of  Sav{‘>')  at  half  maximum 
for  a  square  lattice  with  a  =  0.35  and  SizejrQ  =  20. 

Figure  5:  The  poo  as  function  of  r  for  a  cubic  lattice  of  linear  size  lOOo,  borehole 
radius  5a,  a  =  0.3.  Results  were  averaged  over  50  realizations.  Note  that  when  only 
digable  fragments  are  considered  poo  reaches  1  at  the  maximum  r  the  digable  zone 
reached  in  the  simulations  performed.  However  when  all  fragments  are  considered 
Poo  approaches  1  asimptoticaly  since  the  probability  of  fragment  formation  is  small 
but  not  zero  even  beyond  the  digable  zone. 

Figure  6;  Digable  ‘volume’ (area)  in  two  dimensional  blzisting  on  a  square  lattice 
for  lattice  of  linear  size  1000a,  borehole  radius  50a,  a.  —  0.35.  Results  were  averaged 
over  50  realizations. 

Figure  7:  Digable  volume  in  cylindrical  blasting  on  a  cubic  lattice  for  lattice 
of  linear  size  100a,  borehole  radius  10a,  a  =  0.35.  Results  were  averaged  over  200 
realizations. 

Figure  8;  Digable  volume  in  spherical  blasting  on  a  cubic  lattice  for  lattice  of 
linear  size  100a,  borehole  radius  10a,  a  =  0.4.  Results  were  averaged  over  200 
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realizations. 

Figure  9:  The  excess  digable  fraction  for  the  three  types  of  charges  discussed 
previously.  The  borehole  was  chosen  as  5a  in  all  cases. 
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Abstract 


In  a  series  of  impact  experiments  on  composite 
plates  which  were  aimed  at  testing  the  conjecture 
that  the  perforation  threshold  represented  a  critical 
phenomenon  and  that  Percolation  Theory  could 
eventually  provide  an  adequate  description  of  the 
spallation  or  fragmentation  process  accompanying 
the  perforation  of  anisotropic  layered  materials.  The 
fragments  emerging  from  a  0/90  layup  fiber  epoxy 
plate  target  were  collected  and  photographed.  Image 
analysis  techniques  were  used  to  deduce  the  size 
distribution  of  the  backside  fragments.  It  was  found 
that  the  number  N  of  fragments  having  area  S  may 
be  represented  as  a  power  law  N  ~  1  /  S  ^ .  The 
average  exponent  for  the  experiments  which  were 
performed  at  nine  different  values  of  impact  velocity 
over  the  range  from  300  to  6000  fps  had  the  value  J 
=  2.17  -  18  %  ,  which  lies  between  the  previously 
predicted  theoretical  values  from  Percolation  Theory 
of  Jt20)  =  2.05  and  iTfaD)  =  2.2.  For  the  velocity  V, 
=  397  fps  which  is  close  to  the  perforation  threshold 
velocity,  J  had  the  value  2.03  ,  which  is  surprisingly 
close  to  the  two  dimensional  Critical  Percolation 
Limit. 

introduction 

The  mechanical  properties  of  a  composite  are 
severely  degraded  by  impact  induced  damage.  A 
significant  fraction  of  the  projectile  kinetic  energy  is 
channelled  into  rival  modes  of  damage, ^  i.e. 
deiamination,  fiber  fracture,  matrix  cracking  and 
fragmentation.  Some  disintegration  appears  around 
and  inside  the  penetration  hole,  and  at  the  backside 
of  the  target  spall  is  produced  even  at  velocities  less 
than  the  ballistic  limit.  More  energy  is  consumed 
above  this  limit,  as  intense  fragmentation  emerges 
from  the  target.  This  process  may  be  studied 
quantitatively  by  measuring  the  fragment  size 
distribution. 

Under  a  short  impulse  of  pressure,  a  situation  typical 
of  either  an  impact  or  an  explosion  is  created.  The 
stress  acts  in  a  local  way  and  is  "unable"  to  scan  the 
medium  for  its  weakest  points.  Therefore,  many 
cracks  grow  simultaneously  and  independently  in 
different  regions  of  the  target.  These  conditions  lead 


to  the  well  known  result  (4)  that  the  average  fragment 
size  becomes  smaller  the  higher  the  strain  rate. 
Moreover,  they  provide  the  physical  justification  for 
using  percolation  theory  in  predicting  fragment  size 
distributions.  (2,3) 

Percolation  theory  (1)  has  been  recognized  as  a 
statistical  physics  model  for  critical  phenomena.  By 
applying  it  to  fracture  and  fragmentation  problems  we 
were  able  to  show  that  fragmentation  (as  well  as  the 
coalescence  of  many  microcracks  forming  a 
macroscopic  fracture)  is  a  critical  phenomenon. 
Fragmentation  occurs  whenever  the  density  of  cracks 
reaches  its  critical  value.  That  is  to  say,  that  when 
the  probability,  p,  of  creating  a  microcrack  achieves 
the  critical  value,  p^.  For  crack  densities  higher  than 
Pg,  the  critical  density  of  microcracks,  the  number  of 
small  particles  in  the  population  grow  at  the  expense 
of  the  reduction  of  the  number  of  large  ones. 
Consequently,  the  total  surface  area  of  the  fragments 
is  growing.  It  is  therfore,  expected,  that  the  minimal 
surface  energy  required  for  fragmentation  is  received 
at  Pg  (  below  p^.  a  complete  shattering  of  the  finite 
macroscopic  piece  of  material  is  not  possible). 
Coming  back  to  impact  physics,  the  determination  of 
the  limit  velocity  below  which  perforation  of  a  given 
target  by  a  given  projectile  is  not  possible  continues 
to  be  of  much  interest.  Even  the  definition  of  the 
ballistic  limit  (5)  (the  critical  value  of  striking  velocity 
at  which  the  probability  of  perforation  is  50%) 
suggests  that  stochastic  mechansms  are  at  play. 
Below  the  ballistic  limit,  it  is  expected  that  the  spalled 
fragmented  zone  at  the  backside  of  the  sample  and 
the  shattered  and  penetrated  zone  in  the  front  side 
will  be  well  separated  by  a  zone  of  delaminated 
material  strong  enough  to  prevent  the  projectile  from 
perforating  the  target.  In  the  limit  as  this  intact  zone 
becomes  very  thin,  these  two  zones  will  touch  each 
other,  and  the  sample  will  then  become  disintegrated 
throughout  its  thickness.  It  is  of  interest  to  inquire 
whether  the  critical  condition  for  fragmentation 
supplied  by  percolation  theory  is  able  to  yield 
additional  information  about  the  ballistic  limit  and  the 
physical  conditions  of  its  occurrence. 
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Image  Analysis 

The  population  of  fragments  from  each  experiment 
was  collected  on  the  inner  surface  of  a  cylindrical 
sheet  of  paper  coated  with  a  sticky  adhesive.  After 
the  experiment,  the  cylinder  was  rolled  flat  and 
photographed.  Such  a  procedure  has  to  be  carried 
out  with  caution  in  order  to  avoid  any  optical  noise 
resulting  from  shadows,  unequal  illumination  or 
spurious  particles. 

The  photographs  were  scanned  and  a  digitized  image 
was  obtained.  The  next  step  is  the  division  of  the 
digital  image  into  clusters  of  material  (fragments)  and 
the  background.  This  procedure  is  not  easy  when 
the  gray  levels  of  the  fragments  and  local 
background  overlap.  Also,  it  is  not  easy  to 
distinguish  two  fragments  partially  overlapping  each 
other.  Special  data  treatments  were  adopted  to 
reduce  distortion  and  deterioration  of  the  information 
recorded  in  the  digital  picture,  which  was  scanned  ^y 
a  512  X  512  black  and  white  scanner.  Higher 
resolution  of  the  original  photograph,  of  course,  will 
result  in  higher  final  resolution  and  consequently 
more  data  on  possible  lost  Information  with  regard  to 
the  very  fine  fragments. 

The  measured  fragment  size  histograms  vary 
continuously  in  an  irregular  way.  In  order  to 
overcome  the  effects  of  random  fluctuations,  the 
results  were  accumulated  in  bins  having  unequal 
width.  Each  bin  further  along  on  the  size  axis  was 
made  twice  as  large  as  the  preceding  one.  The  total 
number  N  of  occurences  in  a  bin  was  divided  by  the 
bin  width,  to  obtain  the  average  distribution  density  of 
size  S.  A  linear  least  squares  fit  of  log  N  against  log 
S  was  performed  in  order  to  deduce  the  parameters 
of  the  expected  power  law  relationship.The  results 
are  presented  in  Figure  1,  The  slopes  and  intercepts 
of  the  linear  fits  are  presented  in  Table  1  for  each 
velocity.  ^ 

The  role  of  strain  rate  in  controlling  fragmentation 
was  recognized  by  Grady  in  the  eariy  eighties. 
The  highest  strain-rate  found  during  an  impact  may 
be  estimated  by  the  ratio  of  the  projectile  to  its 
radius,  de/dt  --  Vp  /  Rp  Strain  rates  at  the  backface 
may  be  estimated  by  multiplying  this  value  by  some 
attenuation  function  for  the  particle  velocity.  In  any 
case,  it  is  expected  that  the  higher  the  impact  energy 
is,  the  smaller  will  be  the  average  fragment  size.  It  is 
not  surprising  therefore  that  the  trend  of  intercepts  in 
Table  1  really  represents  the  fraction  of  fine  particles 
in  the  distribution.  The  behavior  of  the  slopes  needs 
further  consideration.  The  mass  of  a  fractal  structure 
having  a  length  L  (L»1)  is  given  by  M(L)  A  L  ^ 
where  D  is  the  fractal  dimension.  At  the  critical  point, 
p  =Pj.  the  value  of  D  is  known.  For  homogeneous 


percolation  in  the  plane,  D  =  1.9  (2D).  In  three 
dimensional  space  D  =  2.5  (3D).  An  additional  result 
of  percolation  theory  is  that  the  number  of 
clusters  having  s  sites  (per  site)  is  given  by  n^  1  / 

,  where  T  =  (d+D)  /  D,  and  where  the  exponent  CT 
has  the  value  S  =  2.05  in  two  dimensions  and  J=  2.2 
in  three  dimensions  (d=3).  All  the  slopes  in  Table  1 
are  around  these  values  with  an  error  not  exceeding 
about  20%.  Furthermore,  in  the  introduction,  a 
possible  relation  between  the  perforation  threshold 
and  critical  fragmentation  was  conjectured.  The 
lowest  velocity  in  Table  1  is  not  far  from  the 
perforation  threshold  found  experimentally  namely 
Vq  =  375  fps.  The  slope  of  the  corresponding 
fragment  ^size  distribution  is  extremely  close  to  the 
theoretical  value  of  J  for  the  two  dimensional  case.  If 
the  fragments  of  the  impacted  composite  plate  have 
a  fixed  thickness  (  as  would  be  the  case  if  it  they 
resulted  from  the  fragmentation  of  previously 
delaminated  plies),  then  the  measured  area 
distributions  would  represent,  in  fact,  two  dimensional 
mass  distributions.  In  a  laminated  material  with 
planar  weaknesses  between  the  plies,  such  a 
hypothesis  appears  plausible.  In  general  isotropic 
solids,  a  true  three  dimensional  behaviour  Is 
expected  with  tT  =  2.2.  It  is  interesting  to  check 
whether  fragment  size  distributions  from  impacted 
brittle  soldis  also  behave  in  the  same  way  near  the 
perforation  limit.  More  work  is  needed  in  order  to 
correlate  other  statistical  Physics  features  either  with 
dynamical  predictions  or  with  experimental  data. 
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Table  l,a:  Number  n(s)  of  fragments  having  area  s  is  a  power 
law:  n(s)  =  A  s""^.  Slopes  and  intercepts  of  the  least  square  best 
line:  Log  n  versus  Log  s,  for  9  impact  velocities  on  (0°/90°) 
(fiber  oriented)  composite.  Slope  extrapolated  to  near  the 
ballistic  limit  equals  to  -2.1. 
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No. 

Velocity 

(ft/sec) 

Slope 

Intercept 

1 

398 

2.28 

316 

2 

1092 

2.15 

49 

3 

1967 

2.45 

792 

4 

2008 

2.40 

580 

5 

3148 

2.69 

861 

6 

4159 

2.63 

817 

7 

4902 

2.87 

1472 

8 

5939 

2.66 

1181 

Table  l,b:  Number  n(s)  of  fragments  having  area  s  is  a  power 
law:  n(s)  =  A  s"'^.  Slopes  and  intercepts  of  the  least  square  best 
line:  Log  n  versus  Log  s,  for  9  impact  velocities  on  (0°/45°) 
(fiber  oriented)  composite.  Slope  extrapolated  to  near  the 
ballistic  limit  equals  to  -2.2. 
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From  steel  ball  impact  on  laminated  composite  (carbon  fiber  reinforced 
epoxy  matrix)  plates  we  have  obtained  the  hole  and  delamination  areas, 
total  ejecta  masses  and  ratio  of  residual-to-initial  velocities  as  functions 
of  Vi-Vc,  the  initial  velocity  (Vi)  less  the  threshold  penetration  velocity 
(Vc ,  the  ballistic  limit).  The  power-law  nature  of  the  dependencies  and 
the  observed  exponents,  both  below  and  above  Vc  ,  indicate  that  plate 
penetration  occurs  as  phase  transition  whose  critical  point  is  Ve  .  We 
found  that  the  energy  absorption  displays  cusp-like  behavior  and 
projectile  velocity  ratio  variation  approaches  universal  power-law 
dependence  typical  of  the  order  parameter.  Confirmation  of  a  critical 
behavior  is  that  the  ejecta  masses  distribute  close  to  the  inverse  2.2 
power  of  the  size,  as  in  a  percolation  process. _ _ 


L  BRIEF  SUMMARY  OF  CRITICAL  PHENOMENA 

Examples  of  critical  continuous  phase  transitions  (Ma  1976,  Callen  1985)  include  the 
phase  change  from  liquid  to  a  vapor  at  the  critical  pressure  and  temperature,  and  the 
vanishing  of  the  magnetization  of  a  ferromagnetic  material  at  zero  magnetic  field  as  its 
temperature  is  raised  past  the  Curie  Temperature.  In  the  case  of  the  former,  it  has 
been  observed  that  large  fluctuations  of  density  develop  which  increase  with  approach 
to  the  critical  temperature.  The  fluctuation  phenomenon  may  be  observed  visually  as 
the  appearance  of  a  critical  opalescence  or  opaque  cloudiness  of  the  fluid  at  the  critical 
point  itself  The  opalescence  is  due  to  the  presence  of  bubbles  of  almost  all  sizes  in  the 
liquid  that  scatter  light  of  all  wavelengths.  This  broad  range  of  sizes  requires  the 
dependence  of  the  distribution  on  a  power  law  of  the  distance  from  the  critical  point 
and,  indeed,  around  the  critical  point,  the  difference  between  the  mean  densities  of 
liquid  and  vapor  have  been  measured  to  vary  as  a  power  law  function  of  the 
temperature  interval  to  the  critical  temperature,  with  a  non-integral  value  for  the 
exponent  of  the  temperature  difference  variable.  The  magnetization  has  a  similar 
behavior,  with  the  quite  astonishing  similarity  that  the  exponents  of  the  power  law 
variations  with  the  temperature  of  these  two  phenomena  are  almost  indistinguishable  to 
within  experimental  accuracy  (Goldenfeld  1992).  The  density  difference  and  the 


magnetization  thus  vary  with  temperature  difference  from  the  critical  temperature  as 
follows:  pf-  pj=  C  (T  -  ,  M  =  C'(T  -  T,)'’-’"*®"'’'. 

The  above  facts  lead  to  the  following  general  statements: 

i)  Fluctuations  of  properties  increase  in  the  vicinity  of  the  critical  point.  The  spatial 
extent  of  the  fluctuations  are  quantified  by  a  reference  length  (correlation  length) 

ii)  A  quantity  called  an  order  parameter  exists  which  is  non-zero  on  one  side  of  the 
critical  point  and  vanishes  on  the  other  side.  The  magnetization  is  such  an  order 
parameter  for  the  Curie  Point  of  ferromagnetic  phenomena.  The  liquid-vapor  density 
difference  plays  the  role  of  order  parameter  for  fluid  systems. 

iii)  Material  properties  obey  power  laws  with  respect  to  the  distance  from  the  critical 
point.  The  specific  heat  at  constant  volume  is  singular  at  the  critical  point  and  displays 
a  cusp-like  variation  on  either  side  of  the  critical  point. 

iv)  The  exponents  in  the  power  law  dependence  of  various  thermodynamic  quantities 
on  the  temperature  difference  to  the  critical  temperature  have  been  conventionally 
given  certain  lower  case  Greek  letters,  as  follows.  The  order 

parameter:  M  oc  (r^  -  ,  The  specific  heat  at  constant  volume:  C,.  oc \T^  -  7j'“ ,  The 

correlation  length:^  k  -  T’l"' .  Experimentally  determined  values  for  these 
exponents  for  various  systems  were  published  in  the  literature  (Goldenfeld  1992). 

v)  The  critical  (universal)  exponents  are  interrelated  by  algebraic  relationships  of  the 
following  typical  forms.  As  a  consequence  of  these  relations,  only  two  of  the  critical 
exponents  are  independent,  the  remaining  ones  being  deducible  from  these  two.  As 
examples  we  cite  the  Rushbrooke  relation:  a  +  2p  +  y  =  2  and  Widom 
relation:  /  =  P{S  -  l)  where  S  is  another  critical  exponent. 

vi)  Mean  Field  Theory  has  not  been  successful  in  exactly  predicting  the  critical 
exponents.  Accepted  modern  computational  approaches  for  computing  the  critical 
exponents  utilize  the  methods  of  the  Percolation  Theory  or  the  Renormalization  Group 
Theory  (Creswick  et  al.  1992).  For  the  Ising  Model  universality  Class  of  Phase 
Transitions,  renormalization  group  techniques  give  the  following  values  for  the  critical 
exponents:/?  =  0.325,  a  =  0.125  while  three  dimensional  Percolation  Theory  for  a 
cubic  lattice  gives:  =  0.4038.  Results  of  various  types  of  theoretical  computations 
for  the  critical  exponents  are  summarized  in  the  literature  (Stauffer  1985). 

n.  DESCRIPTION  OF  THE  EXPERIMENTAL  PROGRAM 
The  experimental  program  consisted  of  firing  one-half  inch  diameter  spherical  steel 
projectiles  at  carbon  fiber  epoxy  matrix  composite  plates  which  were  clamped  between 
steel  plates  containing  a  3  inch  diameter  circular  hole  with  the  design  impact  point 
located  at  its  center.  The  composite  target  plates  were  32  plies  thick.  Stacking 
arrangements  consisting  of  combinations  of  0,  90°  and  ±45°  fiber  orientations  and  also 
in  balanced  spiral  staircase  layups  in  which  the  fiber  direction  increased  in  increments 
of  90,  45,  22.5  or  1 1.25  degrees  from  ply  to  ply.  Ten  inch  diameter  plastic  cylinders 
with  axes  coincident  with  the  shotline  were  lined  with  adhesive  paper  and  placed  flush 
with  the  front  and  rear  target  support  plates  in  order  to  catch  spall  fragments 
emanating  from  the  target  plate.  Circular  end  plates  containing  a  small  central  hole  to 
permit  uninterrupted  passage  of  the  projectile  were  also  covered  with  adhesive  paper 
to  catch  fragments  travelling  with  velocity  vectors  contained  within  a  small  cone  of  the 
axis.  The  population  of  fragments  from  each  experiment  was  photographed  and 
scanned  in  order  to  form  a  digitized  image.  Fragment  size  histograms  were  obtained 
by  using  image  processing  techniques. 


Projectile  velocity  was  measured  ahead  of  the  target  plate  by  means  of  a  set  of  fine  trip 
wires  of  known  separation  distance  and  by  noting  the  instants  of  the  break  of  electrical 
continuity  and  behind  the  target  by  two  sets  of  induction  coils  of  known  location  along 
the  pipe  axis  which  were  wrapped  around  the  outer  cylindrical  surface  of  the  plastic 
pipe.  Thus  the  energy  imparted  by  the  projectile  to  the  target  material  could  be 
computed  (Fig,  1).  Some  experiments  were  performed  without  the  adhesive  paper  and 
the  loose  fragments  collected  in  the  plastic  pipe  were  swept  together  and  weighed  and 
compared  with  the  weight  loss  of  the  target  plates  (Fig.  2).  A  single  set  of  plates 
consisting  of  plates  which  had  been  penetrated  at  different  velocities  were  stained  with 
Gold  Chloride,  pyrolyzed  in  an  air  oven,  deplied  with  a  scalpel,  mounted  on  a  rigid 
backing  and  subjected  to  a  measurement  of  the  areas  of  the  penetration  hole  and  of  the 
envelope  of  the  surrounding  gold-stained  delaminated  area.  This  was  done  with  the 
aid  of  automatic  image  analysis  equipment,  It  was  thus  possible  to  obtain  the  shapes 
and  distribution  of  hole  and  delaminated  areas  on  each  ply  of  panels  penetrated  at 
different  velocities  over  the  range  of  200  to  6000  feet  per  second  (Lair  1991, 
Samad  1994). 

m.  THE  BALLISTIC  LIMIT  AND  CRITICAL  PHASE  TRANSITIONS 
In  ballistic  impact  the  projectile  closing  velocity  rather  than  the  temperature  plays  the 
role  of  marker  for  the  critical  point  of  impact  phenomena.  We  have  therefore  a  ballistic 
limit  velocity  or  penetration  threshold  velocity  rather  than  a  critical  temperature,  and 
closeness  to  the  critical  state  is  measured  as  a  difference  between  the  initial  projectile 
velocity  and  that  at  the  ballistic  limit.  (To  emphasize  the  critical  nature  of  the  ballistic 
limit  we  prefer  the  notion  Vc  over  the  conventional  one  V50)  An  explanation  of  the 
similarity  between  the  initial  projectile  velocity  and  the  temperature  will  be  presented 
further  on  in  this  paper. 

While  critical  points  are  a  phenomenon  of  systems  in  thermodynamic  equilibrium,  most 
ballistic  impact  phenomena  are  considered  to  be  dynamic  phenomena  which  produce 
material  states  that  are  both  irreversible  and  far  from  equilibrium.  Certainly,  gross 
states  of  fracture  and  fragmentation  do  not  appear  to  have  ever  been  observed  to 
reheal  themselves  back  into  an  unflawed  continuum  state.  However,  “reversibility”  is 
not  a  necessary  condition  for  phase  transition  (though  critical  fluctuations  are,  as  a 
rule).  An  example  for  “irreversible”  phase  transition  is  a  percolative  process  which  is 
describable  as  one  with  an  infinitely  large  activation  energy  in  the  reverse  direction. 
More  common  phase  changes,  like  in  ferromagnetism  or  the  Ising  model,  have 
microscopic  (single  site)  activation  energies  that  are  of  the  order  of  the  transition 
temperature. 

Consider  the  trend  of  the  damage  and  the  energy  transfer  ratio,  defined  here  to  be  the 
ratio  of  imparted  to  initial  projectile  energies,  as  the  projectile  velocity  (or  energy)  is 
increased.  At  very  low  velocities,  the  projectile  experiences  a  perfect,  loss-free 
rebound  from  the  panel  and  the  transfer  ratio  remains  at  zero.  As  the  velocity  is 
increased  beyond  the  damage  threshold,  a  crater  is  produced  in  the  panel,  causing  the 
transfer  ratio  to  increase.  With  higher  impact  velocities,  incipient  spall  first  appears 
and  the  transfer  ratio  further  increases.  Below  the  ballistic  limit,  it  is  expected  that  the 
spalled  fragmented  zone  at  the  backside  of  the  sample  and  the  shattered  and  penetrated 
zone  in  the  front  side  (the  crater)  will  be  separated  by  an  intact  or  delaminated  layer  of 
the  material  which  is  strong  enough  to  prevent  the  projectile  from  perforating  the 
target.  In  the  limit  (in  fact,  the  ballistic  limit)  this  intact  layer  becomes  very  thin  and 
the  two  damaged  zones  will  touch  each  other.  That  is  to  say,  the  sample  will  then 
become  disintegrated  throughout  its  thickness. 
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Percolation  theory  has  been  recognized  as  a  statistical  physics  model  for  critical 
phenomena.  By  applying  it  to  fracture  and  fragmentation  problems,  we  were  able  to 
show  that  fragmentation  occurs  whenever  the  density  of  the  microcracks  reaches  its 
critical  value.  That  is  to  say,  when  the  probability  p  of  creating  a  microcrack  achieves 
its  critical  value  pc.  At  and  near  this  point,  the  number  n  of  fragments  (clusters)  having 
size  s  (per  site)  behaves  as  n  cc  1/s^,  where  the  critical  exponent  x  has  the  value  x  = 
2.05  for  a  two  dimensional  isotropic  system  and  x  =  2.2  for  a  three  dimensional 
isotropic  one.  We  have  found  (Murat  and  Anholt  1992)  the  value  x  =  2.1  for  an 
anisotropic  two  dimensional  system  (having  different  crack  probabilities  in  the  x  and  y 
orthogonal  directions).  This  is,  therefore,  the  statistical  physics  predicted  result  for  the 
fragment  size  distribution  at  the  critical  point,  which  is  identified  as  the  ballistic  limit. 
According  to  the  U.S.  Navy  definition  (Zukas  1982)  of  the  ballistic  limit,  the  residual 
velocity  of  the  projectile  is  practically  zero,  thus  producing  a  transfer  ratio  of  unity. 
Any  motion  of  the  target  is  quickly  damped  out  by  the  frictional  resistance  to  motion 
presented  by  the  internal  fracture  (Tamilvanan  1991).  Beyond  the  ballistic  limit,  the 
energy  transfer  ratio  decreased  sharply  because  the  projectile  acquires  an  increasing 
residual  velocity  relative  to  its  initial  value.  At  intermediate  velocities  above  the 
ballistic  limit,  the  transfer  ratio  levels  off  to  a  small  value  which  depends  on  and 
increases  with  the  thickness  of  the  target  panel.  In  the  present  experiments  which  were 
conducted  with  32  ply  thick  plates,  the  ratio  corresponding  to  this  plateau  attained  the 
value  0.2.  At  much  higher  velocities,  the  energy  transferred  divided  by  one  half  the 
projectile  mass  (£„,)  varies  as  E„,  =  CVi"  [V;  =  input  velocity  in  ft/s  C  =  1.33  x  10'^  and 
n  =  1.49  (Altamirano  1990)].  This  causes  the  transfer  ratio  to  decrease  to  zero 
asymptotically  and  inversely  as  the  square  root  of  the  velocity.  The  variation  of  the 
transfer  ratio  around  the  ballistic  limit  (velocity  or  energy)  gives  it  a  cusp-like 
appearance  centered  at  the  ballistic  limit  (Figure  3).  The  first  derivative  of  the  transfer 
ratio  increases  to  a  finite  bound  from  below  at  the  ballistic  limit,  and  experiences  a 
discontinuity  with  a  change  in  sign  there  to  a  finite  negative  bound  from  above,  after 
which  it  increases  rapidly  to  zero.  The  transfer  ratio  derivative  at  the  ballistic  limit 
resembles  the  behavior  of  the  specific  heat  of  a  fluid  at  its  liquid-vapor  critical  point 
IV.  EXPERIMENTAL  SIZE  DISTRIBUTION  OF  COMPOSITE  FRAGMENTS 
Our  experimental  investigation  provides  quantitative  evidence  that  the  ballistic  limit 
represents  a  critical  phenomenon.  This  evidence  comes  from  the  size  distribution  of  the 
ejecta  produced  near  the  ballistic  limit.  The  distribution  density  function  for  fragments 
of  area  s  was  found  to  decrease  with  s  raised  to  a  negative  power.  Results  for  several 
impact  velocities  for  0-90°  fiber  orientation  are  presented  in  Table  1.  The  experimental 
exponents  are  very  close  to  the  values  found  theoretically  for  fragment  sizes  obtained 
from  two  dimensional  bond  percolation  on  a  square  anisotropic  lattice  (x  =  2.1)  or 
percolation  on  a  three  dimensional  lattice  (x  =  2.2)  at  the  percolation  limit  (Jaeger, 
Mayer  etal.  1993). 

We  interpret  our  results  by  supposing  that  the  ejecta  originate  mainly  from  the  target 
region  adjacent  to  the  projectile  at  the  moment  of  passage.  There  and  then,  fracture 
conditions  appropriate  to  the  percolation  limit  obtain,  which  are  reflected  in  the  debris 
size  distribution.  Subsequently  and  elsewhere  in  the  target,  the  crack  density  may 
exceed  the  percolation  limit;  the  resulting  fragments  are  not  ejected  and  are  not 
collected  on  the  adhesive  paper. 

Similar  conditions  apply  to  other  forms  of  fragmentation,  when  the  target  is  not 
confined  and  lie  probably  behind  the  s‘^'^  law  that  appears  to  hold  over  several  decades 
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of  s  in  the  space  debris  distribution  (Englman  1991).  Below  the  ballistic  limit,  the 
fragmenting  materials  of,  and  near,  the  projectile  are  confined  by  the  target  and  the 
fragmentation  occurs  beyond  the  critical  fragmentation  limit  (above  pe).  For  velocities 
above  the  ballistic  limit  (and  p  above  pc)  it  is  expected  to  find  many  more  fine  particles 
(and  less  larger  ones)  in  the  distribution.  These  Percolation  theory  outcomes  agree 
with  the  dynamical  considerations  relating  the  average  fragment  size  with  the  strain- 
rate  (Grady  1982). 

V.  ORDER  PARAMETER  &  CRITICAL  EXPONENTS 

Critical  phenomena  may  be  characterized  by  an  order  parameter,  which'  is  zero  on  one 
side  of  the  critical  point  but  takes  on  a  continuous  range  of  non-zero  values  above  the 
critical  point.  The  order  parameter  either  grows  or  decreases  in  accordance  with  the 
as3Tnptotic  power  law  form  of  the  scaling  hypothesis  near  the  critical  point  with  an 
exponent  p.  For  critical  phenomena  belonging  to  the  Ising  Model  class  P  =  0.325.  For 
percolation  on  a  cubic  lattice  P  =  0.4038.  For  ballistic  impact  and  penetration 
phenomena  various  parameters  at  first  thought  might  be  considered  to  be  candidates 
for  order  parameters,  namely  final  projectile  velocity,  the  mass  of  spall  or  ejecta 
expelled  from  behind  the  target,  and  perhaps  the  volume  of  the  hole  punched  out  by 
the  projectile.  The  degree  of  success  with  which  each  one  of  these  was  found  to 
qualify  by  virtue  of  its  exponent  is  discussed  below.  A  numerical  technique  has  been 
developed  for  simultaneously  extracting  the  limiting  exponent  and  the  ballistic  limit 
from  the  experimental  data  (Mayer,  unpublished). 

(a)  Velocities  Above  the  Ballistic  Limit. 

For  data  which  were  obtained  with  identical,  quality  controlled,  fabricated  panels,  the 
residual  energy  fraction  or  square  of  the  ratio  of  residual  to  initial  the  velocity  was 
found  to  vary  with  distance  to  the  ballistic  limit  raised  to  the  0.649  power.  Thus,  the 
velocity  ratio  will  vary  with  one-half  this  exponent,  or  the  0.3245  power,  which  agrees 
closely  with  the  critical  exponent  0.312  of  the  order  parameter  for  the  Ising  Model 
class  in  the  two  experimental  examples  discussed  earlier  which  also  belong  to  this 
class.  The  best  simultaneous  estimate  of  the  Vc,  the  ballistic  limit  velocity  was  found 
to  be  348  ft/s.  The  evolution  of  ejecta  mass  with  velocity  in  the  form  of  the  scaling 
hypothesis  yielded  a  value  of  351  ft/s  for  the  Vo  and  a  value  of  0.743  for  the  critical 
exponent.  The  ejecta  mass  and  hole  area  depend  on  various  factors  of  the  observed 
energy  and  impact  velocity.  As  shown  in  Table  2,  ejecta  mass  and  hole  area  exponents 
appear  to  be  aberrant  and  display  considerable  variability.  Velocity  ratios  show  much 
less  dispersion  and  approach  the  theoretical  values  for  the  order  parameter  to  a 
satisfactory  degree.  Critical  exponents  for  impact  above  the  ballistic  limit  are  recorded 
in  Table  2. 

(b)  Velocities  Below  the  Ballistic  Limit 

Here,  analysis  of  data  for  rebound  velocity  and  decreasing  part  of  the  delaminated  area 
yielded  exponents  of  0.4467  and  0.8387  respectively.  The  exponent  of  the  velocity 
ratio  appears  to  be  within  12%  of  the  universal  exponent  for  the  order  parameter  as 
determined  from  the  Percolation  theory.  We  cannot  relate  the  exponent  of  the 
delaminated  area  variation  to  any  theory  at  this  time.  The  results  are  presented  in 
Table  2. 

(c)  Analysis 

Denoting  impact  and  residual  quantities  by  1  and  2  (respectively  )  and  the  absorbed 
energy  by  Ea,  and  its  dimensionless  form  by  £’  -E^IE^. 
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Ea  =^1  -^2  =^li] 


By  the  scaling  hypothesis,  and  as  confirmed  experimentally  near  the  V50  =  Vc, 

...-A 


Taking  the  derivative 

dy;-d(vjv)  <  “ 

with  a  denoting  the  exponent  of  the  specific  heat  Cv,  is  equal  to  a  =  -  (2P  -  1)  by  the 
Rushbrooke  relation  a  +  2P  +  y  =  2  provided  y=l,  as  in  mean  field  theory.  Notice  that 
a  yields  the  expected  cusp  like  singularity.  In  the  3D  Ising  model  and  from 
experiments  in  many  materials  a  ~  0.110.  From  our  velocity  ratio  measurements 
a  =  0.1218  (Table  2).  The  fragment  number  data  in  Table  1  provide  a  further  way  to 
obtain  critical  exponents,  through  the  susceptibility  (x)  exponent  y.  %  is  identifiable 
with  the  second  moment  of  s  over  sn(s)  for  all  fragments  up  to  fragment  sizes  of  the 
correlation  length 

Z  =  — TT  fo'"  Since  ;[f  oc  and  g  oc  |r,  -  T\~''  ,where  v  =  0.89  (in  3D),  we 


obtain  from  y  =  (4-t)v  and  x  =  2. 1  (see  Table  1)  the  value  y  =  1.69.  (y  =  1.73  in  3D 
percolation).  The  Widom  scaling  relation  with  5  =  4.8  (from  the  3D  Ising  model)  gives 
P  =  .44  which  agrees  closely  with  several  entries  in  Table  2,  obtained  directly  from 
dynamic  data  (not  from  fragment  distribution),  a  can  now  be  derived  from  the 
Rushbrooke  relation  (without  resorting  to  the  mean  field  theory  value  for  y)  to  find  a  = 
-.44,  similar  to  the  fragment  mass  entry,  though  differing  from  others.  The  exponent 
for  0/45  are  gotten  with  a  similar  analysis. 

VI.  CONCLUSION 

Though  the  exponents  in  Tables  1  and  2,  both  those  directly  observed  and  those 
calculated  from  models,  show  variations,  they  bear  out  the  phase  transition  picture  of 
plate  penetration.  Notably,  the  exponent  P  for  the  velocity  ratio  agrees  in  several  cases 
with  experimental  and  theoretical  values  for  the  order  parameter  in  diverse 
thermodynamic  systems. 

In  a  brittle  material  like  the  epoxy  matrix,  a  natural  cause  for  the  phase  transition  is  the 
coalescence  of  microcracks  coming  from  the  crater  and  spall  regions.  It  has  been 
argued  previously  (Englman  1991,  Anholt  et  al  1993)  that  in  solids  that  are  not 
confined,  at  the  moment  of  fragmentation  the  crack  densities  attain  their  critical  level. 
In  impact  this  happens  at  the  ballistic  limit.  The  identification  of  the  ballistic  limit 
velocity  with  a  critical  point  lays  open  the  field  of  ballistics  to  the  theoretical 
apparatus  of  critical  processes,  including  renormalization  and  finite  size  scaling. 
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Fig  1 .  Energy  Fraction  Absorbed  Versus  Impact  Velocity 
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Fig  2.  Residual  Energy  Fraction  Versus  the  Spalled  Mass. 


Fig  3.  Absorbed  Energy  vs.  Initial  One. 
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MAXIMUM  ENTROPY  PRINCIPLES  IN  FRAGMENTATION  DATA 

ANALYSIS 

R.  Englman 

Soreq  N.R.C.,  Yavne  70600 
Israel 

Abstract 

As  a  first  topic  the  article  treats  observationally  obtained  distributions  of  fi-agment  sizes  and 
shows  that  to  derive  conservation  laws  for  the  fragmentation  process  one  has  to  compute  the 
information  entropy  S  and  then  obtain  the  changes  in  S  when  the  conserved  quantities  are  varied  by 
varying  the  experimental  conditions.  The  criterion  for  the  existence  of  conservation  laws  is  the 
equality  between  the  S-derivative  and  the  Lagrange  multiplier  in  the  Maximum  Entropy  distribution, 
coupled  to  the  non- vanishing  of  the  determinant  formed  out  of  the  second  derivatives  of  S. 
Application  of  the  method  to  ballistically  created  rock  pieces  indicates  conservation  of  the  average  of 
log  (fragment-size)  in  the  break-up  process. 

The  other  topic  is  the  distribution  of  broken  bonds  (as  function  of  their  length  1)  in  a  random 
2D  network,  fragmented  by  pulling.  (The  numerical  simulation  was  performed  by  Dr.  M. 
Lemanska).  The  obtained  distribution  is  in  the  form  of  a  product  of  a  prior  distribution  and  an 
exponential  function  (the  constraint  factor  of  Maximum  Entropy)  whose  exponent  is  linear  in  1. 
This  is  yet  another  case  of  the  constraint  factor  being  a  simple  function  of  the  variable. 


1.  Introduction 


It  is  widely  thought  that  probability  theory,  laws  of  inference  and  predictive  procedures  are 
among  the  most  debated  of  disciplines.  They  have  been  such  since  the  days  of  the  Bemoullis  and 
have  remained  such  up  to  the  present  in  terms  of  fundamentals,  in  spite  of  remarkable  technical 
advances  in  the  use  and  manipulation  of  the  mles  [1-3].  Thus  the  achamed  controversy  between 
ftequentists  and  Bayesians  continues  anabated,  with  each  side  holding  to  its  position  tenaciously. 

The  use  of  Maximum  Entropy  (ME)  as  a  predictive  tool  shares  the  trait  of  being  disputed  in 
principle  and  effective  in  practice  [2].  The  formal  apparatus  of  the  method  is  briefly  set  out  in 
Appendix  1.  Stated  in  verbal  language,  the  method  predicts  that  the  "best"  probability  p(s)  for  the 
occurrence  of  the  value  s  for  a  variable  is  obtained  by  maximizing  the  information  entropy,  subject 
to  constraints  (see  appendices  and  the  next  section  for  quantitative  definitions).  The  method  has 
been  used  in  a  wide  range  of  subjects,  e.g.  statistical  mechanics,  molecular  biology,  pattern 
recognition. 

The  present  author  has  applied  the  method  to  distribution  of  sizes  obtained  by  a  variety  of 
ft-agmentation  processes,  where  the  issue  was:  supposing  you  know  the  constraints  (regularities, 
conservation  laws)  on  the  distribution,  what  p(s)  do  you  anticipate  [4,5]. 

In  this  article  we  address  ourselves  to  two  issues  of  principle  in  the  use  of  the  ME  method: 
First,  given  an  empirical  distribution  p(s)  what  can  you  say  about  the  regularities  (constraints,  laws) 
that  underlie  the  processes  which  have  led  to  the  distribution.  Secondly,  we  treat  the  question  of 
priors  (probabilities  or  weight  that  are  attached  to  each  event  before  receiving  experimental 
information).  Both  subjects  will  be  treated  from  the  viewpoint  of  fragmentation,  though  their 
implications  are  far  wider. 

2.  Deductions  of  and  from  Maximum  Entropy 

In  this  section  we  shall  pose  and  discuss  the  following  two  questions  (application  will  be  given 
in  section  3  and  particular  instances  will  be  worked  out  in  Appendix  2): 

(i)  If  our  fragmentation  experiments  have  resulted  in  a  distribution  of  sizes  s  in  the  form  p(s),  is 
there  a  way  to  establish  that  p(s)  is  a  distribution  with  maximum  entropy? 
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(ii)  What  information  can  be  gleaned  from  a  ME  distribution? 

The  immediate  answer  to  question  (i)  is  in  the  negative:  there  is  no  way  of  learning  from  a 
single,  experimentally  derived  distribution  whether  or  not  it  is  in  accordance  with  Isffi.  For  suppose 

that  experiments  provide  a  set  of  numbers,  e.g.  by  sorting  about  10^  fragments  into  10^’^  bins 
representing  sizes.  Writing 

log  p(s)  =  -g(s)  -  Q  (1) 

(where  Q  is  a  normalizing  constant)  this  has  the  "appearance"  of  a  ME  distribution,  especially  if 
g(s)  has  a  sin:q)le  form,  like  being  linear  in  s,  namely 

g(s)  =Gs.  (2) 

It  then  appears  natural  to  conclude  that  G  is  the  Lagrange  multiplier  and  <s>  is  the  invariant  for  the 
distribution.  However,  this  is  illusory.  To  see  this,  we  note  that  we  can  rewrite  (2)  as 

g(s)  =  Gg  =  GjCs-as^)  +  G2S^  (G2  =  aGj) 

=  G2(s-as^)  +  G2(s^-bs^)  +  G^s^  (G3  =  bG2) 

and  in  many  other  ways  that  would  indicate  other  invariants  and  multipliers  consistent  with  the  data. 
Therefore,  to  summarize,  a  single  set  of  experiments  cannot  establish  the  ME  nature  of  the 
distribution.  To  achieve  this  we  have  to  perform  the  experiments  repeatedly  under  different 
experimental  conditions  that  will  generate  different  distributions.  The  quantities  that  we  wish  to 
vary  through  varying  the  experimental  conditions  are  primarily  the  averages  and  then  the  entropy  S. 
The  analysis  is  as  follows: 

Suppose  that  we  suspect  the  logarithm  of  p(s)  to  be  composed  of  several  linearly  independent 
functions  gjj(s),  gj.(s)...  gj^(s).  (R  cannot  be  larger  than  the  number  of  bins.) 

R 

g(s)  =  ^  Grgr(s)  (3) 

r=l 

We  next  compute  the  entropy  and  the  averages 

S  =  p(s)  log  p(s)  (4) 
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<gr>  =  X  P(s)gr(s) 


from  the  distributions  that  are  obtained  under  several  different  conditions;  in  particular,  under 
conditions  that  the  averages  <gj>  change  infinitesimolly  so  that  one  can  obtain  the  derivatives  of  S 

with  respect  to  <gj>-  For  precision's  sake  we  shall  write  the  derivative  in  the  form  customary  in 

Thermodynamics : 


with  the  averages  in  the  subscript  being  kept  constant. 

We  now  assert  and  shall  presently  prove  that  p(s)  is  a  ME  distributions  if  and  only  if  the 
entropy  derivatives  are  all  equal  to  the  corresponding  coefficients  in  the  exponent,  or 

g;=g,  (7) 

If  these  quantities  are  linearly  independent  functions  of  the  averages,  than  the  averages  are  invariants 
for  the  generation  of  distributions  by  experiment.  If  there  are  some  linear  relations  between  the 
Gj.'s,  as  functions  of  the  means  <gj>,  then  the  invariant  functions  are  superpositions  of  the  gj.’s. 

The  simplest  example  of  this  is  when,  say, 

(G2*=)  G2  =  G4  (=G4  ) 
then  an  invariant  function  is  h(s)  where 

g2(s)  +  g4(s)  =  h(s)  (8) 

If  the  equality  (7)  is  not  achieved  with  the  chosen  set  gj(s)j..gj^(s),  then  one  should  choose 
other  sets  of  g's  to  see  whether  (7)  holds.  Formally,  one  can  expand  the  logarithm  of  p(s)  in  terms 

of  complete  sets,  whose  number  of  components  cannot  exceed  the  number  of  bins  -10^*^,  [6],  but 
this  seems  to  be  prohibitively  laborious.  Clearly,  we  do  not  claim  that  the  observed  distributions  are 
necessarily  ME  distributions. 

However,  if  they  are  such,  then  eq.  (7),  can  be  proved  as  follows: 


Let  the  experimental  conditions  be  described  by  a  set  of  constraints  on  the  functions  gj.(s)  so  that 
their  averages  are  given  and  equal  to  <gj>-  The  ME  distribution  is 

p(s)  =  e-^G,g,(s)  -Q  (9) 

_  g-g(s)  -Q 

having  simplified  by  using  (3).  The  normalizing  constant  is 

Q  =  logXe'"^'^ 

S 

In  (9)  the  coefficients  Gj.  are  so  chosen  as  to  reproduce  the  averages  <gj>.  We  assume  that  this  can 

be  done  uniquely.  (If  not,  we  restrict  our  considerations  to  regions  neighbouring  the  ME  choice  for 
the  Lagrange  multipliers  Gj.).  A  non-ME  probability  p(s)  is  represented  by  the  same  functions  g^(s) 

whose  averages  are  fixed  by  the  values  <gj>  and  an  additional  function  <fg(s).  To  maintain  values 

of  <gj>  as  before  the  coefficients  Gj.  change  to,  say,  G^..  With  adaptation  of  the  notation  in  (3)  we 

write 

g(s)  =]^Grgr(s)  (11) 

r 

and 

=  (12) 

S 

Furthermore 


S  =  X<^r<gi>  +  Q 

r 

s  =  ^  G,,<gr>  +  <5g>  +  Q 


(ME) 

(non— ME) 


(14) 

(15) 


Differentiating  first  the  expression  (14)  for  the  ME  case 


5S 

5<gr> 


•o;=G.+x 


aG 


m 


n  9<gr> 


+  • 


aQ 

a<gr> 


(16) 


However,  from  (10) 
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3<g,>  „  ,  3<8.> 


=  -I 


<gm>- 


3<gx> 


SO  that  the  last  two  terms  in  (16)  cancel  and  eq.  (7)  is  proven  for  a  ME  distribution. 

In  the  entropy  derivatives  5S/5<gj>  for  a  non-ME  distribution  two  terms,  corresponding  to  the 

last  two  terms  in  eq.  (16),  again  cancel  but  one  is  left  with  a  further  term,  as  follows: 


6S  -  3<5g> 

: - =  Gr+T - 


S<gr> 


9<gr> 


The  additional  term  in  (15)  can  be  rewritten  as  the  covariance 

cov(gj.(s)  -  <gj>,  5g,j.(s)  -  <5g,j>)  (20) 

where  the  comma  notation  represents 

5g,j.(s)  =  95g(s)/3<gj>  (2 1 ) 

and  the  differential  operates  on  the  coefficients  Gj^^. 

Though  we  have  not  been  able  to  prove  this  rigorously,  the  cov  are  extremely  unlikely  to 
vanish  for  all  the  indexes  r.  Their  vanishing  would  mean  that  the  coefficients  Gj.  are  unaffected  by 

the  addition  or  the  new  function  5g,  while  still  maintaining  the  averages  <gj>. 

We  have  thus  shown  that  the  entropy  derivatives  are  equal  to  the  corresponding  coefficients  in 
the  probability  exponent  for  a  ME  distribution  and  fail  to  be  such  in  a  non-ME  distribution. 

A  corollary  to  the  proof  answers  question  (ii),  in  that  the  functions  gj.(s)  for  which  (7)  holds 

have  invariant  averages  or,  in  other  words,  reflect  some  physical,  geometrical,  etc.  law  or  constraint 
The  introduction  of  entropy-derivatives  is  a  key  step  in  the  development. 

The  second  derivatives  of  the  entropy 


5<gr>5<gs> 


=  s 


ns 


are  also  important.  It  can  be  shown  that  if  the  determinant  formed  out  of  these  derivatives  vanishes 
identically  as  function  of  the  averages,  then  a  smaller  number  of  the  functions  gj.(s)  can  be  found  (by 

some  linear  superposition  of  the  original  functions)  which  form  invariants.  Thus  the  nonvanishing 
of  the  determinant 

|s.^l*o 

is  a  further  condition. 

It  is  now  possible  to  develop  "Statistical  Dynamical"  formulae  equivalent  to  the  familiar 
thermodynamical  ones.  The  conjugate  variables  in  the  latter  (e.g.,  p  and  -v)  are  replaced  by  the  twin 
set  {Gj,<gj>}  and  the  Legendre  -  transformed  functions  are: 

S({<gr>}),Q({G^})  (22) 

connected  by  the  relation  (14)  as  functions  of  the  averages  [7].  Alternatively,  as  functions  of  the 
Lagrange  multipliers 


Q({GJ)  =  S[<g„>({Gj)]  -  X  <gn>({G,})-G„  (23) 

n 

where,  on  the  right  hand  side,  the  averages  have  been  expressed  as  functions  of  the  multipliers  G^. 
The  relations 


9G, 


=  -<gr> 


complement  the  derivative  relations 


- - =  Gr 

3<gt> 


with  which  we  have  been  busy. 


(24) 


(25) 


Maxwell's  relations  take  the  forms 
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^Gr 

9<gm>  3<gr> 

8<gr>  ^  a<gn^> 

3Gni  3Gr 


(26) 

(27) 


The  equivalent  of  the  Gibbs-Duhem  relations  is  given  by  (17). 

It  has  to  be  emphasized  that  the  validity  of  these  relations  ("Statistical  Dynamics  )  does  not 
depend  on  an  equilibrium  (indeed  the  processes  leading  to  the  distributions  can  take  place  far  from 
equilibrium),  but  only  on  the  continuous  variation  of  the  data  with  the  parameters  of  die  experiment. 
(E.g.,  a  continuous  change  of  the  outcome  of  a  biased-coin  tossing  experiment  with  the  height  of  the 
toss.)  Occasionally  continuity  does  not  obtain:  then  we  have  the  analogy  of  phase  transitions  which, 
however,  are  not  treated  here.  What  is  important  to  note  is  that  we  work  in  a  space  (or  sub-space)  of 
the  data  base,  without  direct  reference  to  the  physical  processes  that  generate  the  data. 


3.  Application  in  Target  Ballistics 


#  In  experiments  by  Matsui  et  al.  [8]  projectiles  of  mild  steel  and  of  rocks  having  velocities 

between  17  and  270  m/s  were  made  to  impact  on  four  kinds  of  rocks  (tuff,  basalt,  granite  and 
dunite)  with  various  shapes  and  sizes.  The  size  distributions  of  fragments  were  plotted  for  twenty 
target  samples  out  of  the  60  or  so  studied.  (We  emphasize  these  numbers  to  focus  on  our  earlier 
assertion  that  a  statistical  analysis  requires  the  availability  of  sets  of  data  under  different  conditions, 
^  rather  than  a  single  set.)  Their  data  are  shown  in  figures  1(a)  -  1(d)  and  a  numerical  representation  is 

exhibited  in  Table  1.  The  authors  of  Ref.  [8]  fit  their  data  to  an  inverse  power-law  dependence  of 
distribution  p  on  the  size  s  relative  to  the  total  target  size.  In  our  terminology  this  has  the  form 

^  p(s)  =  s’^  exp(-Glogs  -  Q)  (28) 

where  Q  is  a  normalizer,  s'^  is  the  prior  and  G  the  Lagrange  multiplier.  In  the  choice  of  the  prior  we 
^  have  followed  Jeffrey's  argumentation  [9],  based  on  the  scale  invariance  of  the  distribution.  The 

multipliers  G  (=  a  - 1  in  the  notation  of  [8])  are  given  in  table  1;  they  are  not  constant  but  depend 
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# 


# 


# 


mildly  Gogarithmically)  on  the  kinetic  energy  of  the  projectile.  For  thoroughly  fragmented  targets  G 
ranges  between  1  and  2.  (The  values  of  G  are  close  to  the  fractal  exponents  of  Turcotte  [10]  found 
for  a  wide  assortment  of  fragmentation  processes).  As  shown  in  the  previous  section,  the 
verification  of  (28)  being  a  ME  distribution  is  contingent  to  satisfying 


3<logs> 

To  form  derivatives  from  numerical  data  we  require  data  in  which  two  values  of  <logs>  are 
close.  Inspection  of  Table  1  shows  that  there  are  four  pairs  of  targets  in  which  <logs>  fall  within  .03 
of  each  other.  The  following  table  (Table  2)  shows  that  for  three  pairs  (two  in  tuff  and  one  in 
granite)  equation  (29)  is  satisfied  to  within  ten  per  cent,  while  in  the  remaining  pair  (in  basalt)  the 
discrepancy  is  large  (even  the  signs  of  the  two  members  differ). 

One  might  thus  conclude  that  only  in  tuff  and  granite  is  the  distribution  that  of  ME,  but  not  in 
basalt.  In  fact,  however,  the  discrepancy  in  the  last  is  due  t  o  another  cause:  the  minimum  sizes  (the 
lower  cut-off)  in  B 13  and  B6  differ  (column  4  in  Table  2). 

Now,  the  cut-off  s^  can  be  incorporated  in  the  ME  formalism  by  a  constraint  of  the  form 

I  H(Sq-s)  p(s)ds  =  0 


where  H  is  the  Heaviside  step  function.  (The  corresponding  Lagrange  multiplier  is  infinite).  The 
differentiation  in  (29)  implies  [as  in  eq.  (6)]  the  constancy  of  all  other  invariants  (or  invariance  of 
Sq).  It  appears  that  violation  of  this  requirement  has  led  to  the  failure  of  satisfying  equation  (29). 

In  summary,  then,  the  experimental  data  of  Mitsui  et  al.  are  representable,  by  eq.  (28),  as  ME 
distributions  with  the  invariant  being  <  logs>.  It  is  not  clear  at  this  stage  what  it  is  in  the  physical 
process  which  brings  about  fragmentation  that  causes  <logs>  to  be  invariant.  However,  only  very 
recently  has  a  suggestion  been  made  concerning  a  physical  cause  for  self-similar  distribution  [11] 
and  apparently  much  work  needs  to  be  done  before  relations  between  distributions  and  their 


governing  processes  will  be  established. 
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4.  Priors  and  Simple  Constraints 

A  source  of  difficulty  in  the  use  of  ME  are  the  priors  n(s)  (described  formally  in  Appendix  1). 

They  are  weights  ascribed  a  priori  (that  is  to  say,  without  regard  to  observation)  to  each  event  and 
are  analogous  to  "measure"  in  set  theory.  In  the  words  of  Skilling  [2]:  "The  prior  is  the  Lebesgue 
measure  associated  with  [the  variable],  which  must  be  given  before  an  integral  can  be  defined".  An 
alternative,  probabilistic  interpretation  identifies  priors  with  subjective  beliefs  as  regard  to  outcome 
of  experiments  or  with  relative  reliance  attached  to  different  events.  Only  in  quantum  mechanics  is 
the  choice  of  priors  ('equal  weight  to  all  states’)  incontrovertible,  and  even  there  only  for  discrete 
Hilbert  spaces;  whereas  for  distributions  in  a  classical  or  in  a  continuous  domain  prior  choices  are 
flexible.  Atempts  to  find  "best"  priors  were  based  on  invariance  principles  [9]  or  the  weighting  of 
several  expert  opinions  in  view  of  new  evidence  [12].  [Bayesians  tend  to  use  "priors"  in  the  sense  of 
hypotheses  that  are  already  there  and  need  not  to  be  invented,  but  have  to  be  verified  by  observation 
[13]). 

A  different  understanding  of  priors  and  indeed  of  the  ME  method  was  provided  by  Levine  and 
coworkers,  mainly  in  the  realm  of  molecular  dynamics  [14],  and  extended  recently  to  classical 
processes  (fragmentation)  [15].  In  this  approach,  the  factorization  of  the  distribution  function  P(s) 
into  two  factors,  as  shown  in  Eq.  (A8)  of  Appendix  1,  is  based  on  the  philosophical  principle  that 
Nature  is  generally  kind  enough  to  make  the  second  exponential  factor  short  and  simple.  This 
happens,  e.g.,  in  a  canonical  distribution,  where  the  prior  is  the  multiplicity  or  statistical  weight  of 
the  energy  shell  and  the  exponent  is  simply.  -  (Energy/kg  x  temperature).  Molecular  collision  and 

photofragmentation  experiments  also  lead  to  simple  exponential  factors  (called  also  constraint  or 
likelihood  factors)  [14],  provided  that  the  statistical  or  probabilistic  information  is  condensed  in  the 
prior.  Then  the  constraint  factor  incorporates  the  characteristics  of  the  process  which  led  to  the 
distribution.  This  process  is  almost  invariably  very  complex  and  it  is  only  through  Nature's  so  far 
unfathomed  benevolence  that  the  experimental  data  can  be  approximated  by  constraint  factors  that 


are  simple. 
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Staple  exponendal  facto,  rule  many  empirical  disuiburions.  e.g.,  gauasians  and  log-nonnal 
They  alao  work  for  disnibudona  .ha.  arise  from  numerical  staulariona.  E.g.  diauibudoua  of 
„„  networks  can  be  well  repreaenKd  by  d.e  product  of  a  probabiUatic  prior  and  a 


polygons  in  cellular  networks  can  be 
simple  constraint  factor  [16,17].  In 


the  following  section  we  employ  a  ME  distribution  function, 
mstraint  factor,  to  describe  the  statistics  of  broken  bonds  in  a 


consisting  of  a  prior  and  a  simple  constraint  factor,  to  aesuiu.  u*. 
simulated  Voronoi  -  network  that  has  been  fragmented  by  the  appli^^^^^ 


5.  Fragmentation  in  a  Disordered  Network 

To  model  fraoure  ta  granular  mMcrials  with  geomemcal  disorder  a  pyo  dimensional  model  of 
roUds  consisdng  of  polygonal  grains  has  been  constac«d  by  Voronoi-«sselladon  [18],  This  means 
mat  (UP  m  2<«0)  poin.  were  chosen  randomly  in  a  plane  and.  by  joining  up  the  perpendicular 
Pisecors  ofbonds  Uta.  link  neighbouring  poin.,  dte  plane  was  divided  into  polygons  of  drfferen. 
sires  and  shapes  (Figure  2).  Ihe  average  number  of  sides  of  dte  polygons  is  6.  A  dteory  of  stde 
dismibudon  in  2D  and  in  3D  has  been  given  in  119],  The  present  Voronoi-tessellabon  ts  a 
reasonable  approxtaadon  m  disordered  muldphase  stactures  containing  irregular  grams.  More 
pr^isely,  it  is  turned  superior  to  a.—  models  in  which  disord«  enmrs  dmough  randomness 
of  bond  snengths  in  a  geomettically  regular  (e.g.  square  or  triangular)  ladice  [20^1]. 

In  this  simuladon  fracture  is  brought  about  by  applying  either  a  fixed  stress  or  a  fixed 
elongadon  (constant  sbain)  at  two  paraliel  edges  and  calculadng  the  elongadon  in  each  bond  for  a 

harmonic,  central  force  model  of  springs. 

one  or  mom  bonds  break  when  their  calculated  reladve  extensions  (defined  as  dte  rado  of  the 

extensiontothe  Origin^  bond  lengthlexceedacertain  value.  Then  the  stass  or  strain.^ 

is  mised  undl  newer  bonds  bmaks  (and  in  fact  they  do  so  with  ^ater  facility  since  dte  bmken  bonds 
no  longer  support  the  stmss).  The  procedum  is  repeated  undl  one  gets  a  condnuous  domain  of 
broken  bonds  cutting  across  the  sample,  which  is  then  considered  fractured. 

Several  parameters  of  the  ftacture  process  are  computed,  such  as: 

(a)  dte  leas,  sness  or  stain  which  causes  fracmm  (the  fiacuue  stmss  or  stain). 

(b)  dtenumberofbondsdtatambrokenateachsmgeofdtefracutm  process. 
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(c)  the  location  of  the  broken  bonds  as  function  of  their  length, 

(d)  the  distribution  of  broken  bonds  as  function  of  their  length  (Fig.  3)  under  condition  that  a  fixed 
elongation  AL^  (in  the  x-direction,  say)  was  imposed  on  two  edges. 

A  simple  ME  analysis  of  the  data  shown  in  Fig.  3  is  as  follows: 

The  distribution  q(5l^)  of  elongations  51^^  of  bonds  in  the  x-direction  (II  ALj,)  can  be  obtained 
from  the  conservation  condition  for  the  elongation  component  51^^^  of  the  i'th  bond: 

AL,  =  M-'X5i;  (31) 

i 

where  M  is  the  number  of  instances  that  the  elongation  has  been  applied  The  ME  distribution  then 


follows: 


N  being  the  total  number  of  bonds. 

The  fracture  criterion  used  in  the  simulation  was  for  the  elongation  51^  of  the  i-bond  of  original 


length  Ij 


•^>ef(=0.2)  (33) 

Assuming  that  the  bond-directions  between  points  in  a  Voronoi-tessellation  are  isotropic,  we  can 
replace  the  above  relation  by  the  following  relation  (true  in  the  mean)  for  the  elongations  in  the  x- 
dircction 


5ljj 

Then  the  number  of  fractured  bonds  of  length  1  is 
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Nf(l)  =  Nod)  J  dxg(x) 

kf  /V2 


(34) 


=  Npd)  exp[-  lef/V2  (MAL^^/N] 


The  factor  Nq(1)  is  the  number  of  bonds  of  length  1  created  originally  and  is  shown  in  Fig.  3  as 
the  full-lined  histogram.  serves  as  the  prior  for  the  distribution  of  broken 

bonds.  The  percent  of  broken  bonds  [  Np/N^  x  100]  is  shown  in  Figure  4  for  three  simulations  in 
which  elongations  AL^^  (~  0.04)  were  applied  to  the  end  points  (numbering  about  20)  several  times 

(about  5  times)  until  the  bond  structure,  containing  about  800  bonds  (N  =  800),  became  fragmented. 
The  least  square  fits  of  the  data  to  exponential  curves  (seen  in  Fig.  4)  were  reasonably  accurate.  In 
terms  of  our  earlier  discussion  of  the  roles  of  the  prior  and  of  the  constraint  (or  likelihood)  factor  in  a 
ME  distribution  the  present  simulations  confirm  the  expectation  that  the  exponent  in  the  constraint 
factor  takes  a  simple  form  (linear  in  the  bond  lengths  1),  in  spite  of  the  complexity  of  the 
firagmentation  process. 

The  fitted  exponents  in  Fig.  4  are  about  10  1.  Comparing  this  with  the  exponent  in  Eq.  (34), 
using  also  the  appropriate  values  of  the  symbols  (M,N,  latter  predicts  a  much  faster 

decay  (30 1)  of  the  distribution  curve. 

6.  Conclusion 

Fragmentation  is  a  complex  phenomenon  whose  outcome,  the  fragment  distribution,  is 
expected  to  reflect  this  complexity.  The  ME  approach  has  been  used  in  the  past  to  impose 
regularities  on  the  distribution  [4,5,14].  In  this  article  we  set  out  to  uncover  regularities  firom  the 
data  .  The  key  to  this  effort  are  the  derivatives  of  the  entropy  with  respect  to  the  invariants.  The 
procedure  leads  to  a  "Thermodynamics"  of  the  fi:agment  -  (or  any  other)  data.  Application  has  been 
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made,  to  relatively  simple  case  of  ballistic  fragmentation  of  rock  targets,  where  only  one  invariant 
<logs>  is  in  evidence.  For  future  uses  of  the  procedure,  when  more  than  one  invariant  is  expected, 
further  development  of  the  method  is  needed  to  make  it  practical.  This  appears  possible. 

A  prior  -  constraint  function  decomposition  of  distributions  has  been  earned  out  in  this  article 
for  the  broken  bonds  in  a  disordered  network  at  the  moment  of  its  falling  to  pieces.  The  aim  was  to 
see  whether  the  constraint  function  is  simple,  as  anticipated  in  the  ME  treatment.  This  was  found. 
Additional  interst  is  attached  to  possible  changes  in  the  distribution  as  the  fragmentation  threshold  is 
approached. 
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Appendix  1. 


ThR  formalism 


In  deriving  the  distribution  p(s)  of  fragment  sizes  s  the  ME  method  maximizes  the 
information  entropy  (or  missing  information)  given  by 

S  =  -  ^  p(s)  logp(s)  (Al) 

S 

One  allows  for  R  observational  inputs  or  constraints  in  the  form  of  averages  whose  values  are  given: 

<gr>  =  ^  gr(s)p(s)  (r=l,...,R)  (A2) 

S 

Then  one  maximizes  with  respect  to  p(s)  the  quantity 

S-XGr<gr>.  ^^3) 

r 

where  Gj.  are  Lagrange  multipliers,  whose  role  is  to  ensure  that  with  the  solution  p(s)  the  constraints 
are  satisfied.  This  solution  is 


p(s)  =  exp  -^Grgr(s)-Q 


where  Q  normalizes  p.  Exp  Q  is  also  the  partition  function. 

When  there  are  "priors"  n(s),  discussed  in  section  4  that  are  not  uniform,  the  previous 
quantities  change  in  the  following  way; 


S  =  ^  n(s)p(s)log  p(s) 

s 

<gi>  =  ^  n(s)p(s)gf(s) 


and  the  distribution  is 


P(s)  =  n(s)p(s) 


=  n(s)exp  -  ^  Gj.  g/s)  -Q 

r 
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Appendix  2 


These  illustrate  different  cases  that  are  treated  in  the  text  in  a  general  manner.  (For  convenience 
averages  are  shown  with  bars  over  the  variable). 

Suppose  that  we  find  that  the  experiments  yield  a  distribution  p  of  a  quantity  z  that  can  be 
represented  as 


p(z)=e‘^'""°2" 

Several  possibilities  arise. 

(1)  Two  constraints  (in  addition  to  normalization). 


These  are  on  the  mean  and  the  variance,  in  the  form 


gi  =  z,  gi  =  a 


g2  =  z^  g2=«^ 


Then  p  in  (A9)  has  the  form; 


p(z)  = 


sp  -  ■j(z-a)^/((7^-a^) 


(AlO) 


From  this  we  find  immediately 


Gi  =  -a(a^-a^)  \  G2  =  'j(a^-a^)  ^ 


Q  =  laW-a^)+|log[aW)] 


(All) 


S  =  GjZ  +  G2Z  +  Q 

=  |  +  |log[(27t(o^-a^)] 


The  entropy  is 
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The  partial  derivatives  of  S  with  respect  to  the  averages  are: 


=  G*  =  -  djicr  -  z) 


9S  _  *  1  ,  2  2n-1 

_  =  G,=  jto-a) 


which  are  respectively  equal  to  Gj  and  G'2;  This  verifies  that  (AlO)  is  a  ME  distribution. 


(2)  One  constraint 
(2a)  g'  =  z  +  z2,g'  =  a  +  (r^ 

but  no  constraint  on  z  and  2?  separately. 

In  (A9)  we  have  now  G2^=  G2. 

The  ME  distribution  satisfying  the  one  constraint  is 


p(z)  27tfo^  +  a  +  — 


fix)m  which  the  coefficient  in  the  exponent  is 


(2  i 

G'=  2  a^  +  a+-7 
V  4 


S  =  ^  +  —  log  2rt^o^  +  4 


^p[-(z+i 


5) 


Differentiating  this 


^  =  G'  =  [2(o^.a.i)f  =  G' 


(A12) 


so  that  (A  12)  is  also  a  ME  distribution. 

However,  (AlO)  satisfying  two  constraints  has  a  lower  entropy  than  (A12)  satisfying  one 
constraint  (since,  in  general,  adding  constraints  makes  the  distribution  less  uniform).  Indeed 


AS  1,2  (=  Sq.  -  Sgj.g^)  =  \  log  +  a  + 


(A13) 


This  entropy  difference  is  plotted  in  Fig.  5  as  function  of  o  for  different  choices  of  cP"  and  is  seen  to 

be  non-negative  (for  all  values  of  a  that  make  sense:  P  <  (P  ).  Note,  however,  that  the  entropy 
difference  can  be  negative  for  certain  choices  of  the  priors,  in  particular  of  non-uniform  ones) 

(2b)  g'  =  z\  g=P 

without  constraining  z. 

p(z)  =  (27taV^^  exp^-  j  =  exp(-  gP  -  Q)  (A  14) 

G'  =  (2a^)"\  Q  =  ^  log(27ia^) 

S  =  -t-  log(27ca^) 

^.G*  =  (2aV'  =  G 
Bg' 

demonstrating  the  ME  nature  of  the  distribution  (A  14).  The  entropy  difference  compared  to  the  2- 
constraint  distribution  is  again  positive: 

AS  =  log{oHP-P)  >  0 

(3)  We  now  consider  a  non-ME  (NME)  distribution  and  show  inter  alia  that  partial  derivatives  of 
the  entropy  with  respect  to  averages  are  not,  in  general,  equal  to  the  corresponding  coefficient  in  the 

exponent  of  the  distribution  p,  i.e.G-  Gj 

(3a)  Assume  that  there  is  just  one  constraint 
g  =  z,  g  =  a 

and  no  constraint  on  P.  Suppose  further  that  the  observation  can  be  represented  (somewhat  to  our 

surprise)  by  the  distribution 

/  \  — Gz— Cz^  —  Q 

p(z)  =  e 

To  satisfy  normalization  and  z  =  a,  one  must  have 
p  =  (2jt/C)”^^  exp[-C(z-a)^] 


(A15) 
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whence 


G  =  2Ca 


Q  =  Ca^  +  -log(2)i/C) 


S  =  --log(23r/C) 
r 

°  -  2C  aq 

oz 


(A16) 


Since  the  variation  of  C  with  a  or  z  is  not  prescribed  G  will  not  in  general  equal  G.  Should 
one  find  experimentally  such  dependence  of  C  on  a  that  G  =  G  is  satisfied,  then  a  new  constraint  on 

z^  has  been  discovered. 

(3b)  Constraint  on  g  =  z^:  g  =  o^,  no  constraint  on  z. 

Suppose  that  observations  yield  the  non-ME  distribution 


-Az-Gz  “Q 

p  =  e 


Satisfaction  of  the  constraint  leads  to  the  quadratic  in  X  -  G”^, 


S  =  ^  +  jlog\2n(o^-jAX^ 


BA  d<j2'"aA 


The  last  non-equality  is  again  due  to  the  arbitrary  dependence  of  A  on  the  mean  i 
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Table  1.  Analysis  of  rock  fragmentation  data  [8].  Computed  averages  are  over  observed 
distribution  points  z^.  z  =  s/s^qj^2  variable  normalized  to  total  mass.  A  prior  ^ 

was  used  to  calculate  the  cross-entropy  (=£  Pn  log  Pi/^)-  z  is  the  minimum  normalized  size  that 

n  ^ 

was  sampled. 


Sample 

<:z>xlO^ 

<log  z> 

Cross-entropy 

ZQ(minimum 

Q 

size)  X  10 

• 

Tuff  T5 

1.642 

-6.524 

-.280 

1.270 

T19 

2.610 

-6.069 

-.246 

1.953 

T13 

2.533 

-6.095 

-.283 

1.953 

• 

T14 

4.245 

-5.622 

-.226 

3.008 

T4 

1.760 

-6.410 

-.328 

1.465 

T16 

2.474 

-6.098 

-.278 

1.914 

• 

T6 

6.518 

-5.086 

-.362 

5.625 

T7 

6.820 

-5.055 

-.307 

5.625 

• 

Basalt  B16 

4.961 

-5.452 

-.307 

3.750 

B13 

4.416 

-5.650 

-.1964 

2.851 

B3 

2.767 

-5.992 

-.298 

2.187 

• 

B5 

4.341 

-5.583 

-.297 

3.281 

B6 

4.004 

-5.627 

-.285 

3.125 

Granite  15 

4.172 

-5.607 

-.222 

3.047 

• 

D 

4.172 

-5.386 

-.219 

3.789 

11 

3.927 

-5.633 

-.255 

3.047 

K1 

3.293 

-5.853 

-.200 

2.344 

• 

Dunite  D4 

2.564 

-6.081 

-.263 

1.953 

D3 

2.302 

-6.171 

-.272 

1.797 

D2 

3.536 

-5.704 

-.324 

2.969 
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Table  2.  Entropy  derivatives  computed  from  data  points  [8]  for  pairs  of  rocks,  compared  with 
corresponding  coefficients  G[T  =  mff,  B  =  basalt,  IG  =  granite].  Also,  smallest  relative  fragment 
size  Zq(z  =  s/SjQjgjjX  =  <  logz>) 


Rock  1/Rock  2 

dS/dx 

G 

Zo(xlO^) 

T19/T13 

1.42 

1.7/1.4 

1.95/1.95 

T6/T7 

1.94 

1.94/1.84 

5.615.6 

B13/B6 

-4 

0.87/1.28 

2.85/3.1 

IG5/IG3 

1.25 

1.4/1.5 

3.05/3.05 
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Figure  Captions 


Fig.  1.  Fragment  size  distributions  for  four  target  types  after  high  velocity  impact  obtained  by 
Matsui  et  al,  [8] .  In  the  log-log  plots  N  is  the  number  of  fragments  per  unit  volume  size 
s.  Sj  is  the  target  volume.  The  different  symbols  relate  to  different  samples,  whose 

data  appear  in  Table  1. 


Fig.  2.  Voronoi-tessellationofaplane.  (Part  of  the  plane  is  shown).  The  lines  are  perpendicular 
bisectors  of  the  bonds  (not  shown),  joining  neighbouring  points  (thrown  at  random). 
The  broken  lines  are  bisectors  of  those  bonds  that  have  been  fractured  by  the  application 
of  tensile  strain  at  the  boundary  points. 


Fig.  3.  Histogram  of  broken  bonds  as  function  of  the  original  bond  length  at  fracture.  The 

represent  averaging  over  several  runs.  The  criterion  for  bond  rupture  is  that  the  relative 
elongation  51/1  of  bonds  exceeds  0.2.  The  distribution  N^G)  of  original  lengths  is  also 
shown  by  broken  lines  and  the  percentage  of  bonds  broken  Nj,(1)/NqG)  by  dotted  lines. 
In  our  model  short  bonds  have  a  preference  for  being  broken. 


Fig.  4.  Fraction  of  broken  bonds  (  percent  of  the  original  number)  as  function  of  original  bond 
length  at  the  moment  when  the  network  becomes  fragmented.  Data  points  from  three 
simulations  are  shown  together  with  fitted  exponentials. 


Fig.  5.  The  entropy  difference  in  (A13).  It  is  non-negative  since  the  distribution  in  (A12) 
has  less  constraint  (more  disordered)  than  that  in  (AlO). 

The  three  curves  are  for  =  1  (broken  line),  5  (full  line),  1  (dotted  line). 


Bond  number 
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BROKEN  BONDS  (percent) 


100.00 
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PART  6 


Transport  Treatment  of  Crack  Population  In  Finite  Medium 
Part  II  :  Damage  Behavior  in  Different  Zones. 

Miemanska,  R.  Englman  and  Z.  Jaeger 
Soreq  Nuclear  Research  Center,  Yavne  81800,  Israel 


Abstract 


Extensive  graphical  results  are  given  of  a  non-linear  transport  equation  approach  (derived  and 
described  in  part  I  of  this  work)  for  crack  population  in  a  spherical  medium  surrounding  a  charged 
borehole.  The  central  quantity  numerically  computed  by  time  integration  is  the  "damage"  or  the  total 
volume  associated  with  the  cracks.  Typically,  intensive  cracking  (large  damage)  appears  near  the 
borehole  and  close  to  the  surface  (spalling),  where  the  detonation  pulse  becomes  (upon  reflection) 
tensile.  The  extents  and  quantities  of  damage  in  these  two  zones  are  given  as  functions  of  the  sample 
size,  initial  crack  distribution,  attenuation  length  and  the  form  of  the  pulse. 

An  interesting  result  arising  from  this  non-linear  method  is  the  disappearance  of  the  spall  zone  in 
spite  of  the  relatively  strong  tensile  part  of  the  wave.  This  phenomenon  is  caused  by  the  increased 
attenuation  of  the  tensile  pulse  (on  its  way  to  the  surface)  due  to  the  large  number  of  cracks  (either 
initially  present  or  formed  later). 

1.  Introduction 

This  work  is  a  continuation  of  the  work  called  here  Pan  I,  in  which  we  have  presented  the 

motivation  and  aims  of  our  research.  The  transport  equation  of  interest  has  been  described 
there.  The  notations  and  definitions  hold  also  here.  In  order  to  facilitate  the  reading  of  this  paper  we 
recall  some  of  the  basic  equations  in  SecJ2. 

The  numerical  simulation  is  performed  for  spherical  geometry  as  in  Part  I,  but  the  model 
considered  here  is  more  complicated.  In  this  work  a  spherical  detonation  source  of  small  radius  a^  is 

assumed,  which  is  concentric  with  the  sample  sphere  of  radius  R^.  The  outgoing  pressure  wave 

consists  of  compressive  (c,  positive)  and  tensile  (t,  negative)  rectangular  parts  with  variable 
amplitudes  and  lengths.  The  wave  is  moving  from  the  center  to  the  boundary  and  is  reflected 
inwards  to  the  center.  The  reflection  turns  the  tensile  and  compressive  parts  into  each  other. 

Both  Part  I  and  Part  n  of  this  work  should  be  regarded  as  an  attempt  to  apply  a  statistical- 
deterministic  method,  the  transport  equation,  to  the  treatment  of  damage  problems  in  a  spherical 
medium.  To  test  the  usefulness  of  this  method,  we  investigate  here  the  dependence  of  the  final 
damage  on  (i)  the  initial  crack  distribution,  (ii)  the  sample  size,  (iii)  the  pressure  wave  and  (iv)  the 
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wave  attenuation  factor.  We  also  examine  the  effect  of  the  sample  dimension  on  the  extent  of  the 
cmshed  regions  (  around  the  borehole  and  near  the  surface). 

2.  Assumptions  and  Definitions 

We  work  with  a  pressure  wave  of  the  form: 


Fig.  1 .  Depiction  of  the  pressure  wave.  The  symbols  have  the  following  meanings: 

:  strength  of  compressive  and  tensile  pressure  respectively. 

Wq  w^:  width  of  compressive  and  tensile  parts  of  the  pressure  pulse. 

Next,  we  recall  some  basic  equations  described  in  detail  in  Part  I.  Let  r  be  the  crack-  tip  position, 
£  its  long  dimension,  v^  is  the  crack  velocity  postulated  to  be  isotropic  and  constant .  The  non-linear 
transport  equation  for  crack-tips  motion  at  time  t  is: 


3T  -  f  r 

^  +  Vf  n  gradT  -  P(r,At)  -  K*(£'^£,  r,t,T,T’)  d^'  -i-  K"  (£^£',  r,t,T,T')  d^'  =  0  (1) 

^  J('<i  Jc>t 


where  T  (r,^,t)  dr  d^  is  the  number  of  crack-tips  in  the  ranges  (r+dr,r),  {£+^,£)  at  time  t . 

The  source  term  P: 

P(r,4t)=7i  I  a(r,t)  -  G;  1  (GoTp^)"'  (2) 

where  £^  is  the  mean  crack  length,  the  amplitude  of  detonation  pressure  and  tp  stress-pulse 
time.  They  are  assumed  to  be  unity,  is  crack  activation-rate  coefficient  (i=c,t),aj  is  a  threshold 
pressure  for  crack  growth,  both  are  dependent  on  the  stress  type. 


The  scattering  terms  are  similar  to  the  source  term.  Explicitly 
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- K*  lU'  + K-  «■  =  Xi  I  c(r,t)  -  Oi  I  (o<,y,r'  -f  T(r/,t)s 


A.T  d/-] 


where  -  R -a„  is  the  maximum  length  of  cracks,  X-  is  crack  growth-rate  coefficient  (i=c,t)  and 

V/  w  •  1 


A=j^”“d£j^  !  4^^“ 


T(r,^',t)  df 


-^nux/4 

Jo 


T(r,i’,t)  d£' 


The  pressure  function  aCr.t) 


a  (r,t)  =  a*  F(r,t) 


where  a*  =  for  a  compression  and  tensile  pressure  respectively.  The  function  F 

determines  the  attenuation  process  in  term  of  the  path  integral  of  the  local  damage  D(r,t). 


F  (r,t)  =  -^  exp  -  f  D(i’,t)  dr’(Ro  5  ^ 

^  Ja 


where  5  is  a  medium  dependent  attenuation  factor.  The  damage  at  point  r  and  time  t  is  defined  by 


D  (r,t)  =  r  T(r,At)  d£ 
Jo 


The  final  damage  is  also  defined  by  Eq.(7)  with  T(r,^,t|^^),  where  t^^^^  is  the  travelling  time  of 
the  pulse  from  r=aQ  to  the  boundary  and  backward. 

The  initial  crack  distribution  is  assumed  to  be  uniform  in  the  medium  and  is  given  as  a  function 
of  crack  length  £  : 

T„(r,4  =  T(r,f,0)=-^e-^^^o  (8) 

This  presents  the  initial  condition  for  our  transport  simulation. The  parameter  a  may  be  varied, 
but  in  general  in  this  work  it  is  put  a=1.0. 


The  initial  damage  is  denoted  by  D^.  Assuming  the  above  uniform  crack  distribution,  is 
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constant .  From  Eqs.(7)  and  (8),  it  is  obtained 


(70 


The  ultimate  mean  damage  is  defined  as  the  ratio  of  total  damage  to  the  sample  volume, 

V 

sample* 

^mean  =  ^  D(r)  dr  / 

o 

where  D  (r)  is  the  final  damage  at  the  point  r  and  time  t 

max 

We  next  characterize  the  spall  as  the  damage  near  the  external  boundary  by  its  depth  and 
Its  peak  value.  Similarly,  for  the  damage  zone  around  the  borehole.  (For  illustration  we  show  Fig.2) 


Fig.  2 

For  the  smdy  of  items  (i)-{iv)  in  the  previous  section,  the  wave  parameter  values  are  chosen 
after  some  preUminary  computations,  with  a  view  of  yielding  results  that  compare  well  with 
observations,  as  follows: 

a  =  10.0a  w  =0.1  w 

^  U  c  p 

a^  =  0.01  Gp  Wj  =  0.9  Wp 

where  a^  is  assumed  to  be  unity,  (see  Table  1  of  Part  I). 

The  parameters  in  the  above  expressions  were  given  the  values: 


^Q-1.0,  Xpl.O,  Yc=l-0  for  compressive  wave  and  Yt  =75.0  for  a  tensile  one. 
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3.  Effect  of  The  Initial  Crack  Distrihution  on  The  Final  Damage. 

The  influence  of  the  initial  crack  distribution  on  the  final  damage  is  examined  for  two  values  of 
the  sample  radius:  =  5.0^^  and  =  8.5-f^j.  Troughout  this  work,  the  radius  of  the  spherical 

charge  is  kept  constant  3^=  0.5/^^.  The  pressure  pulse  is  also  the  same,  with  parameter  values  as  in 
Sec.2. 

The  initial  crack  distribution  is  changed  through  variation  of  the  parameter  a,  see  Eq.  (8). 

We  used  the  following  values:  a=0.01,  0.1,  0.5,  1.0,  2.0,  5.0  and  10.0.  Consequently  the  initial 
damage  density  varies  from  0.03944  to  39.44  for  RQ=5.0fjj  and  from  0.0574  to  57.4  for 

Ro=8.5fo  (See  Eg.(7')  ).The  range-of  a  values  is  sufficiently  large  for  discussing  the  question  of 
interest. 

Due  to  the  similarity  of  the  computational  outcomes  in  both  cases,  only  the  ones  for  R  =  5M  are 

0  0 

presented  in  Figs.  3,4,5.  For  better  estimation  of  the  obtained  results,  the  differential  damage 
produced  by  the  explosion,  which  is  defined  by  D  (r,t)-DQ,  is  shown  in  Fig.  4  and  the  relative  one 
D  (r,t)/D^  in  Fig.  5. 

The  following  conclusions  can  be  drawn.The  ratio  D(r,t)  /D^  grows  as  the  value  of  a  decreases. 

In  other  words,  a  smaller  initial  damage  results  in  getting  later  a  relatively  larger  one.  This  is  clear, 
because  the  presence  of  cracks  causes  the  attenuation  of  the  outgoing  pulse.  Indeed,  it  is  found,  that 

for  sufficiently  large  a  (^),  for  both  R^  /  £^=  5  and  8.5  cases,  the  spall  zone  does  not  exist. 

For  sufficiently  small  a  (^.1),  the  peak  and  the  depth  of  the  spall  are  independent  of  a.  The 
same  result  is  observed  for  the  inner  damage  zone.  See  Fig.3.  That  is  clear,  because  the  initial  value 

of  the  function  F(r),  [see  Eqs.(6)  and  (7)  ]  is  small  and  has  no  influence  on  the  computational 

results.  In  the  limit,  as  a  -4o,  the  peak  of  the  damage  near  the  center  reaches  a  rnim'mal  value  and  the 
peak  of  the  spall  attains  a  maximal  one. 
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•  Effect  of  the  Sample  Size  on  the  Final  Damfl^re  nistrihiifinn 

In  this  section  the  radial  behaviour  of  the  final  damage  for  various  sample  sizes  will  be 
^  considered.  To  investigate  this  quantity,  is  varied  between  l.Of^  and  10.5^^.  Other  parameters 

are  kept  fixed  and  have  the  values  as  in  Eqs.(lO)  and  (11).  The  parameter  a  is  also  fixed  and 
assumed  a=1.0.  The  results  are  presented  in  Figs.  6,7,8. 

First,  we  consider  the  damage  caused  by  the  outward  pulse,  as  shown  in  Fig.6.  The  outward 
pulse  near  the  boundary  is  weak,  and  its  contribution  to  the  damage  is  small.  In  this  figure  the  spall 
does  not  appear  because  it  is  created  only  by  a  strong  enough  reflected  (inward)  pulse. 

^  The  damage  peak  near  the  charge  increases  as  grows,  but  levels  off  for  sufficiently  large  R^,  to 

become  independent  of  the  sample  radius.  This  may  be  explained  by  Eq.  (7*). 

Writing  for  simplicity,  oc,/=l  in  (7’)  and  noting  that  the  maximum  crack  length  t  =R -a  we 
_  obtain: 
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=  f  e~^  d£  =6  —  (  ^ ^  e  ^  df  6  for  Rg  — >  «> 

•'o  •'V^ 


(7") 


That  means,  for  sufficiently  large  R^,  the  value  of  the  function  F(r),  Eq.(6),  is  independent  of  the 
sample  size. 

The  final  damage  plotted  in  Figs.7  and  8  is  caused  both  by  the  outward  and  returning  pulses.  For 
,  the  width  and  the  peak  of  the  inner  damage  zone  rise  quickly.  An  interplay  of  two 

factors  influences  this  phenomenon.  One  is  the  behavior  of  the  initial  damage  explained  above. 
Generally  speaking,  the  second  one  is  the  ratio  of  the  pulse  length  to  the  sample  size.  For 
R^=1.0  the  two  lengths  are  equal  .The  outward  and  the  reflected  pulses  act  simultaneously  at  all 


points  r,  and  the  contribution  of  the  returning  wave  is  small,  because,  in  fact,  only  is 

effective.  Note,  that  these  a  differ  in  sign.  If  the  pulse  length  is  smaller  than  the  sample  radius,  the 
two  pulses  act  simultaneously  only  in  a  limited  region  of  points,  r,  near  the  boundary  and  the 
contribution  of  the  returning  pulse  is  important 

For  sufficiently  large  R^,  R^>4.5f^,  the  peak  and  width  of  the  inner  damage  zone  are 

independent  of  the  sample  size.  This  may  be  explained  by  an  interplay  between  two  factors.  The  one 

is  the  behavior  of  the  outward  pulse  considered  above  and  the  second  one  is  that,  the  reflected  pulse 

is  weak  for  large  R  and  therefore  its  contribution  is  small. 
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Fig.  6.  Radial  behaviour  of  the  damage  produced  by  the  outward  pulse  for  various  sample  sizes. 
For  R^=8.5,9.5,10.5  only  the  damage  region  around  the  borehole  is  depicted. 
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Fig.  7.  The  final  damage  distribution  as  a  function  of  the  sample  size. 
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Fig.  8.  The  final  damage  distribution  as  a  function  of  the  sample  size. 

The  behavior,  of  the  peak  and  width  of  the  spall  is  different  from  the  above.  The  peak  spall,  D^, 
dependence  on  the  sample  size  R  may  be  approximated  as  a  parabola  D  =4.73R  ^+41.61R  - 
39.76,  in  the  range  2.5^q5.5.  It  rises  quickly,  when  R^^  increases  from  1.0/^  to  about  4.5^^,  and 
attains  its  maximal  value  for  R^  ~  4.4^^^.  For  R^  >4.4^^.  the  peak  decreases  and  for  sufficiently 
large  R^  ,  the  spall  does  not  exist.  This  is  caused  by  the  pulse  behavior,  which  becomes  attenuated 
during  its  passage  due  to  material  damage. 

In  most  cases  considered  above,  two  damage  zones  are  observed.  One  exists  near  the  borehole 
and  the  other  near  the  boundary.  Between  these  two  zones  there  is  a  region  of  low,  almost  constant 
damage.  When  is  this  value  of  damage  equal  to  the  initial  one?  What  are  the  relations  between  the 
extent  of  ±e  region  of  constant  damage,  the  inner  damage  zone  and  the  sample  size  Rq  ?  In  order  to 

answer  these  questions,  we  let  R^  vary  between  1.5.^^  to  13.0^^.  The  parameters  and  a  are 


given  by  Eqs.(lO)  and  (11).  The  value  of  the  parameter  5,  see  Eg.(6)  is  chosen  to  assure  the 
existence  of  these  three  seperate  zones  as  S=4.225.  The  results  of  many  computations  are  collected 
in  Fig.9.  For  0<Rq^.25^q  ,  the  damage  at  any  point  r  is  greater  than  the  initial  one,  curve  (a),  For 

4.25^g  <Rq<5.0^q,  two  seperated  damage  zones  (internal  and  spalled)  and  a  region  of  initial  damage 

are  observed.  For  this  value  of  5  and  when  two  damage  zones  are  obtained:  an  innpir 

damage  zone  and  an  initial  one.  (The  spall  disappears).  According  to  the  above  results,  the  curve  (c) 
indicates,  that  the  width  of  the  inner  damage  zone  is  indifferent  to  the  values  for  sufficiently  large 

samples.  Of  course  the  ranges  of  cited  above  are  dependent  on  the  value  of  5.  Similar  provisos 
may  be  made  in  connection  with  the  effects  of  other  parameters,  as,  for  example,  the  pulse  ones. 


Fig.9.  The  extent  of  the  damage  zones  versus  the  sample  size.  Curve  (a):  the  damage  at  any  point  is 
greater  ±an  the  initial  one,  D^,  the  extent  of  the  damage  zone  is  equal  to  R^-Curve  (b):  the  extent  of 

region  versus  R^.  Curve  (c):  the  width  of  the  inner  damage  zone  near  the  borehole  versus  R^. 
Curve  (d):  the  width  of  the  spaU  versus  R^. 


5.  Width  and  Volume  of  Spalled  Regions  for  Samples  of  Various  Sizp 

During  the  study  of  this  problem  the  sample  radius,  R^,  is  varied  and  other  parameters  are  fixed 

as  in  Eqs  (10)  and  (11).  The  spall  zone  consists  of  regions  of  various  damage  value.  Their  widths 
for  various  R^  are  depicted  on  Fig.lO.  In  Fig.ll  the  spall  zone  volume  and  in  Fig.l2  the  mean 

damage,  versus  R„  are  drawn. 

®  mean’  o 
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Fig.  10.  Width  of  the  subregions  of  the  spall  zone,  having  damage  greater  than  D.(D.= 
20,25,27,28,29,30,35)  versus  the  sample  size  R^.  For  example:  for  Ro=4.5  the  width  of  the 
spaU  subregion  of  D>35  is  0.52  the  width  of  subregion  of  D>25  is  1.5  In  the  case  of 
Ro=8.0 there  does  not  exist  a  spall  subregion  of  damage  greater  than  35. 


Fig  1 1  The  volume  of  the  spaU  zone  subregions  versus  R^,  for  D.=25,30,35. 


In  this  section  the  pressure  wave  remains  rectangular  and  consists  of  a  compresion  and  tensile 
parts,  but  the  value  of  the  parameters  a  and  T  are  varied  in  order  to  examine  their  influence  on  the 
final  damage.  The  sample  size  is  Ro=5.0^q  and  ^^=0.51^.  Three  cases  will  be  considered  below. 

a)  The  values  of  the  parameters  a^,  Tp  and  the  length  w^.,  w^  of  the  compression  and  tensile,  parts 
of  pulse  remain  as  in  Eqs.(lO),  (1 1).  The  amplitude  of  the  tensile  wave  varies,  takes  the  values: 

-0.001g^,-0.01Gjj,-0. IcJ^  j,  and  -  I.Og^.  It  is  noteworthy  that  the  role  of  the  tensile  wave  in  the 
damage  process  is  more  significant  than  that  of  compressive  one.  Moreover  in  the  returning  phase 
of  the  wave  the  roles  of  G^  and  are  changed,  g^  -^g^.  and  cyg->Gj.  In  the  light  of  these  remarks,  the 
results  plotted  on  Fig.  13  are  easily  understood.  The  large  peak  of  the  damage  near  the  borehole,  in 
the  case  of  G^  =  I-Og^,  is  caused  by  the  outward  pulse.  During  the  passage  to  the  boundary  this 
pulse  becomes  so  weak  (due  to  the  created  large  damage),  that  the  reflected  one  is  unable  to  form  a 
spall.  For  considerably  smaller  value  of  G^,  ±e  effect  of  the  outward  pulse  is  weak  and  the 
contribution  of  the  returning  one  is  more  important 
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b)  In  the  second  case  the  values  of  the  w  ,  w_,  w,  remain  as  in  (a).  But  here,  a,  is  fixed  and  its 

P  C  I  t 

value  is  a^=0.01aQ.  varies  and  takes  the  values:  0^=1. 5.0a^,  lO-Oa^  and  20.0aQ.  The 
obtained  results  are  shown  in  Fig  14. 


c)  In  this  case  the  parameters  a^,  /  Wp  and  /  Wp  are  kept  fixed  and  their  values  are  as  in 

Eqs  (10)  and  (11).  The  parameter  Xp,  the  pulse  time,  varies.  Note  that  since  the  pulse  velocity  is 
constant,  the  pulse  length  varies  accordingly  with  Xp  values,  Wp=0.5  £^,  1.0^^^  and  2.0^^^.  For 
larges  Xp  both  x^  and  x  ^  are  large,  which  explains  the  results  given  in  Fig.  15. 


Fig.  15.  Final  damage  for  various  wave  length. 
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From  the  results  of  the  cases  (a)  and  (b)  it  may  be  concluded,  that  the  ratio  plays  an 
important  role  in  the  spall  -  formation.  As  was  observed  in  our  computations  for  spalling,  this  ratio 
should  be  sufficiently  large,  cr^/  <7^  >100.  It  is  clear,  that  the  values  which  cause  spalling 

depend  on  other  parameters,  especially  on  (crack  growth)  and  6  (attenuation).  On  the  basis  of  our 
computational  experience,  we  found  that  even  though  <5^  /  a^<100  a  large  enough  y^  creates  a  spall. 
Same  questions  remain  still  open,  such  as  whether  there  exist  a  value  -  range  of  generating  a 

spall  ?  Whether  there  exists  an  optimal  value  of  <7^  ^  ^^^h  the  spall  peak  attainc  a 

maximum?  In  order  to  find  answers,  additional  work  should  be  performed.  Similarly,  it  is 
anticipated,  that  there  exists  an  optimal  pulse  length,  for  which  the  spall  peak  achieves  a  TnaYimnm 

7.  Effect  of  the  Attenuation  Factor  Value  on  the  Final  Damapp 

The  attenuation  constant  5,  which  appears  in  Eq  (6),  plays  a  remarkable  role  in  decreasing  of 

the  pressure,  a(r,t),  during  the  pulse  passage  over  the  sample.  The  factor  5,  characterizes  a  material 
property.  It  is  important  to  examine  its  effect  on  the  final  damage. 

In  this  section  all  parameters,  except  5,  have  the  values  given  by  Eqs.(lO)  and  (11). 

First,  we  examine  the  effect  of  various  5  values  on  the  radial  distribution  of  the  final  damage.  The 
computations  were  performed  for  two  samples  having  =  3.0 and  =5.0^^^.  For  both  samples 

^  ^  Siven  the  values:  75.0,  175.0,  250.0,  1000.0,  10000.0. 

In  the  limiting  case,  5  F(r)-^l/r.  See  Eq.  (6).  For  5M^.  5  is  given  the  values:  2.0, 

10.0,  50.0,  200.0  and  2000.0.  Furthermore,  the  limit  case,  F(r)=  l/r,  is  also  considered.  The 
results  obtained  for  both  two  cases  ^  similar  and  therefore  only  the  one,  for,  R^=3.0^ij  are  depicted 
in  Fig.  16. 

The  damage  near  the  sample  center  rises,  as  the  5  value  grows.  For  sufficiently  large  5,  the 
damage  decreases  quickly  as  a  function  of  the  position  and  ±e  spall  does  not  exist.  To  understand 

these  results,  we  note  that  the  function  F(r)  increases  with  5  (See  Eq.  (6)).  Therefore  for  larger  5, 
the  pulse  is  only  weakly  attenuated  near  the  borehole. 
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Fig.  16.  Effect  of  6  value  on  final  damage. 

Now  we  shall  examine,  what  are  the  relations  beween  the  inner  damage  width,  the  spall  width, 

the  extent  of  the  initial  damage  region  and  the  value  of  the  5  parameter.  5  will  be  varied  and  all  other 

parameters  will  be  kept  fixed  at  the  values  given  by  Egs.(lO),  (11).  Three  sample  sizes  were 

considered,  =5.0£^,  and  R^=8.5.^^  .  The  case  of  R^=5.0^^  will  only  be  discussed  in 

detail,  because  the  results  obtained  for  other  values  of  R  are  similar 

o 

To  estimate  the  influence  of  the  initial  damage,  the  parameter  a  takes  the  values  2.0, 1.0  and  0.1. 

The  results  for  a=1.0  are  ploted  in  Fig. 17,  where  three  curves,  are  drawn.  Curve  (a)  gives  extents 
of  initial  damage  regions.  On  the  curve  (a')  are  shown  the  cases,  for  which  a  spall  does  not  exist.  On 
curve  (a")  the  case  with  spall.  Curve  (b)  gives  extents  of  the  inner  damage  zones  and  curve  (c)  the 

spall  widths.  For  example,  for  5=4.5  the  extent  of  region  is  0.45  ^^(A),  the  width  of  spall  is 

1.25  (B)  and  of  inner  damage  is  2.8  ^  (C).  The  sura  of  these  extents  is  4.5  £  .  Note-  R  -a  =4  5 

For  0<5<4.75  the  inner  damage  zone  and  the  region  of  the  initial  damage  are  observed.  The  spall 
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exists  only  for  3.65<5<4.75.  For  6=4.75  the  extent  of  the  initial  damage  region  is  equal  to  zero.  The 
inner  damage  zone  and  spall  zone  are  coupled  and  their  widths  achieve  a  maximum.  The  damage 
distribution  for  this  situation  may  be  sketched  as  in  Fig.  1 8. 

For  5>4.75  the  damage  at  any  point  r  of  the  sample  is  greater  than  the  initial  one,  and  its  behavior 
was  discussed  at  the  beginning  of  this  section.  Another  interesting  point  is  the  limit,  5— >0.  Then 
F(r)^0,  see  Eq.(6),  and  the  pulse  is  inactive  . 


Fig.  17.  Extent  of  damage  zones  as  a  function  of  5  value  for  a=l  .0. 


1.0  2.0  3.0  4.0  5.0  6.0 

Fig.20.  Extent  of  damage  zones  as  a  function  of  5  value 
curves(a),  (b)  and  (c)  is  as  in  Fig.  17. 
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At  last  we  should  note,  that  the  curves  on  Fig.  17  are  approximately  straight  lines.  This  fact  is 
remarkable,  because  our  transport  equation  (see  Sec.2)  is  not  linear  and  the  5  parameter  is  included 
in  an  exponential  term. 

On  Figs.  1 9  and  20  the  results  for  a=0.1  and  a=2.0  are  given.  They  indicate,  that  ±e  value  of  5, 
for  which  a  spall  exists,  increase,  as  a  grows.  Similar  results  are  obtained  for  and  R^= 

8.5^0- 

Additional  results  shown  in  Fig.21  are  worthy  of  note.  Each  curve  in  this  figure  presents  a 
scaling  relation  between  those  values  of  a  and  6,  that  causes  a  spall  of  given  width.  The  curves  on 
Fig.21  refer  to  the  spall  width  equal  to  0.2^^,  0.5^^,  sample  size  R^=5.0^^  and 

a  =0.5-^  .  The  relation  is  linear.  Similar  results  are  obtained  for  R  =3.0^^. 

0  0  0  0 


Fig.  21 .  Relation  between  a  and  5  for  given  spall  widths,  which  are  denoted  on  the  curves. 
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8.  Conclusions 


The  great  number  of  graphical  data  in  this  article  arises  from  the  extensive  computational  work 
earned  out  throughout  several  years  (on  an  IBM  machine).  Occasional  repetitiveness  is  excused  by 
the  need  to  show  that  the  results  are  systematic  and  what  are  the  tendencies  in  the  results. 

The  transport  equation  approach  to  crack  population  growth  is  only  one  of  the  possible 
treatments  of  a  complicated  stochastic  problem:  to  us  it  appears  neither  less  nor  more  justified  than 
any  other,  and  for  a  space  and  time-dependent  description  it  may  be  uniquely  workable.  There  are, 
admittedly,  a  great  number  of  free  parameters  in  the  approach,  whose  values  are  not  easily  fixed  in 
any  concrete  situation.  Our  aim  was  to  elucidate  the  wide  range  of  phenomena  that  can  occur  with 
some  (not  unphysical)  choice  of  the  parameters.  Yet  more  parameters  can  be  varied  (as  the  crack 
velocity  or  introducing  inhomogenity  in  the  initial  crack  distribution). 

Although  the  resultes  presented  in  this  paper  are  related  to  the  samples  defined  in  Secs.  3-7,  it  is 
possible  to  reach  some  conclusions,  which  seem  to  present  the  general  behavior  of  damage  systems. 

(a)  In  a  given  sample,  a  pressure  wave  creates  a  spall  zone,  if  the  ratio  odot  is  larger  than  some 
threshold  value. 

(b)  For  a  given  pressure  wave  in  a  given  medium,  the  damage  peak  in  the  spall  zone  increases 
for  Ro  <  R’^'o  and  decreases  for  Ro  ^  R*o.  For  sufficiently  large  size  of  the  sample  Ro,  the 
spall  does  not  exist.(  In  the  case  of  Sec. 4,  R*o=4.4) 

(c)  For  a  given  pressure  wave  in  a  given  medium,  the  damage  peak  and  the  width  of  the  inner 
damage  zone  increases  with  the  sample  radius  Ro.  Both  level  off  and  became  independent  of 
the  sample  size  for  sufficiently  large  Ro.  (See  Sec.4) 

(d)  For  a  given  pressure  wave  and  sample  size  Ro,  if  a  >  a*,  where  a*  is  a  sufficiently  small 
number,  the  damage  peak  and  the  width  of  the  inner  damage  zone  and  the  spall  region  are 

independent  of  the  a  value.  In  the  limit,  as  a— »0,  'the  peak  of  the  damage  near  the  center 
reaches  a  minimal  value,  and  the  peak  of  the  spall  attains  a  maximal  one.  (  See  Sec.3 ) 

(e)  For  a  given  pressure  wave  and  sample  size  Ro,  if  ot  >  ct*,  where  cc*  is  sufficiently  large 
number,  the  spall  zone  does  not  exist.  (  See  Sec.3  ) 

(f)  A  pressure  wave  and  a  sample  of  radius  Ro  are  given.  For  5  >  5*,  is  a  sufficiently  large 
number,  the  spall  zone  does  not  exist  and  the  damage  peak  of  the  inner  damage  zone 

increases  with  the  5  value.  In  the  limit,  5 — >co,  the  damage  achieves  a  maximal  value. 
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